はじめに
有限要素法に代表される非構造メッシュベースの数値解析手法は複雑な任意形状を解析できる点で、CAE (Computer Aided Engineering) の中核的な技術として、構造解析、流体解析、熱伝導解析などに広く用いられます。計算精度や計算時間の観点に立つと、二次元解析や構造要素(シェル要素)を用いた解析では四角形メッシュ、三次元解析では六面体メッシュが適用されます。三角形メッシュや四面体メッシュは、Delaunay 分割により自動メッシュ生成の手法が確立されています (例えば Shewchuk et al., 2016)。一方、四角形メッシュ、六面体メッシュを精度よく自動生成する手法は確立されておらず、これらメッシュ分割を行う際には手作業による修正を含む多くの労力を要しているのが現状です。
科学計算総合研究所(Research Institute of Computational Science Co. Ltd., RICOS)では、四角形メッシュ、六面体メッシュを精度よく自動生成する技術について研究開発を進めています。四角形メッシュ・六面体メッシュを高精度に生成する手法として、frame field に基づくメッシュ生成が期待されています (Ray et al., 2006, K$\ddot{\textrm{a}}$lberer et al., 2007, Bommes et al., 2012, 2013b)。この手法は対象領域があらかじめ三角形メッシュ・四面体メッシュでメッシュ分割できることを前提として、
- 四角形メッシュ生成では $\pi/2$ の回転対称性
- 六面体メッシュ生成では立方体対称性
を有する frame の姿勢分布を最適化することで、四角形メッシュ・六面体メッシュをそれぞれ抽出する方法です (Bommes et al., 2009, Ray et al., 2016)。
この記事では、frame field によって二次元四角形メッシュを生成する技術について紹介します。
Frame field による二次元メッシュ生成
Frame の最適化
Frame field によるメッシュ生成手法は、対象領域に対して $\pi/2$ の回転対称性を有する frame の角度分布を最適化することで、四角形メッシュを抽出する方法です。Frame は言うなれば、四角形メッシュの「もと」のような役割があります。
はじめに、frame の角度分布の最適化結果を見てみましょう。図 1 に、扇形領域に分布する frame field の例を示します。この例では、外形線に沿った frame の角度を境界条件として与え、内部の姿勢を最適化しています。
図 1:扇形領域に分布する frame field の例。
この最適化結果から、特異な構造があるかどうかを観察します。赤線は、frame field の流線です。この赤線をたどってみると、対象領域のほとんどは 4 つの流線で囲まれていることがわかります。一方中央にひとつだけ、4 つの流線で囲まれていない領域(三角形の領域)が見つかります。この領域を特異(singularity)といいます。
また、ある点周りの frame field が特異である場合、その点を特異点(singularity point)といいます。任意の形状で完全四角形分割を考えるとき、どうしてもある点から伸びる辺の数が 4 以外になってしまうことがあると思います。実際にそのような点が特異点と対応しています。詳細な定義は例えば Vaxman らの論文 (Vaxman et al., 2016) などに譲ります。
図 1 に示す最適化された frame の分布を見つめると、局所的に四角形メッシュの構造(四角形メッシュの「もと」)が隠れているのがわかります。つまり frame field による二次元メッシュ生成とは、この frame 分布の情報を使って四角形メッシュを抽出する方法です。
Integer Grid Maps による四角形メッシュ抽出
最適化された frame 分布の情報から四角形メッシュを抽出するために、integer grid maps を利用します。Integer grid maps の概略は以下の通りです。
- Integer grid とは座標 ($u_1$, $u_2$) を考えたとき、$u_1$, $u_2$ どちらかの値が整数値をもつもの
- 対象領域の境界(外形線)を integer grid に乗せて、格子の情報を対象領域に引き戻せば、四角形メッシュが抽出できる(図 2)
- 任意の対象領域から四角形メッシュを抽出するには、対象領域に切込みが必要(図 3)
- 切込みの終点は特異点でかつ、座標 ($u_1$, $u_2$) のどちらも整数値でなければならない(図 3)
- 切込みを入れたことで生じる境界のペアが integer grid 上で綺麗に一致しなければならない
- 切込みは特異点または境界上の点を結ぶ曲線で、自己交差せず、他の切込みとも交差しない
対象領域の境界(外形線)を integer grid に乗せるような変換ができたとき、integer grid の情報を元の領域に逆変換することで四角形分割が得られます。例えるならば、integer grid 上で得られる格子のテクスチャを元の領域に綺麗に張ることができれば、その格子を四角形分割とみなすことができるということです。
一方このような変換は、境界の座標が必ず整数値になるような制約と、切込みが integer grid 上で綺麗に一致する制約の2つを満たすことが必要です。前者の制約は、上記で説明した通りです。後者の制約は、切込みを入れたことで生じる境界のペアが、integer grid 上で整数値の並進移動と 90度の整数倍の回転で一致することを指しています。
そこで、上記のような四角形分割を与える変換を得ることができれば万事解決です。しかしながら、この変換を得るのは容易くありません。なぜならば、前述した整数変数に関する 2 つの制約を与えられているため、解くべき問題が整数混合最適化問題になるからです。整数混合最適化問題は、連続変数と整数変数が混在する最適化問題で、NP 困難であることが知られています。
提案手法
整数制約条件の緩和
整数混合最適化問題を貪欲的に解く場合、整数制約条件の数を $M$ として、$O(M)$ 回の線形方程式の求める回が必要になります (Bommes et al., 2009, 2010)。また、メッシュ生成時に得られる整数混合最適化問題を貪欲的に解く場合、整数制約を緩和することが影響してメッシュ抽出できない解を得ることが指摘されています (Bommes et al., 2013a)。
そこで、整数制約を何らかの形(ここでは線形制約)に緩和して、整数制約の数を削減できないか考察してみることにします。メッシュ生成に必要な制約は、境界の座標が必ず整数値になるような制約と、切込みが integer grid 上で綺麗に一致する制約の 2 つでした。前者の制約は完全四角形メッシュ生成のために重要な制約ですのでこれ以上の緩和は難しそうです。後者の制約はどうでしょうか。
後者の制約は、integer grid から引き戻してきたときに定まる四角形メッシュのエッジを予想して切込みを入れると、integer grid 上で境界が全て $u_1$軸と$u_2$軸に一致します。具体例を図 4 に示します。後者の制約のお気持ちは、切込みによって生じた境界が integer grid 上で綺麗に一致してほしいということでした。今考えている切込みは既に$u_1$軸と$u_2$軸に一致しているので、後は切込みによって生じた境界の長さをinteger grid 上で一致させる制約を与えれば、後者の制約は満足することになります。
また、この切込みは三角形要素のリメッシュによって表現することとしました。従来手法では三角形要素の中心をノードとするグラフ構造でなされていた frame の最適化ですが、リメッシュのため、有限要素離散化によって frame の最適化を行いました。このとき当該の制約は$u_1$と$u_2$の線形制約になります。具体的な切込みの与え方や線形制約の与え方は、著者らの論文 (森田ら, 2020) に拠ってください。
このように、切込みに一工夫すると integer grid 上で整数値の並進移動と 90度の整数倍の回転に関する整数制約の個数を削減することができます。整数制約条件の数を削減すると、メッシュ生成に必要な計算時間が削減できることが予想されます。次章ではその有効性を検討することにします。
図 4:Integer grid から引き戻してきたときに定まる四角形メッシュのエッジを予想して切込みをいれると、integer grid 上で外形線がu1軸、u2軸に沿う例。
数値例
整数制約条件の数を削減すると、メッシュ生成に必要な計算時間が削減できることが予想されます。そこで、提案手法によって得られた四角形メッシュの有効性を検証するために、4 つの基本形状に対してメッシュ品質と計算時間を比較します。メッシュ精度の指標として、aspect ratio と scaled Jacobian を用います。Aspect ratio は 1 以上の値をとり正方形では最小値 1 をとる指標です。Scaled Jacobian は、-1 から 1 の値をとり正方形では最大値 1 をとる指標です。
既存手法との比較:メッシュ品質
三角形形状、円形状、五角形形状、リング形状の 4 つの領域に対し、既存手法と提案手法を比較します。本数値例では既存手法として、frame field と integer grid を計算するライブラリ libigl と、integer grid から 四角形メッシュを抽出するライブラリ QEx を利用しました。
図 5 に、それぞれのメッシュ分割例を示します。表 1 に、メッシュ分割例の aspect ratio の最大値、平均値、最小値を示します。提案手法と既存手法を比べると、提案手法は円形以外の例で aspect ratio の最大値を低減できており、品質の良いメッシュ分割が得られたことがわかります。表 2 に、メッシュ分割例の scaled Jacobian の最大値、平均値、最小値を示します。提案手法と既存手法を比べると、提案手法は scaled Jacobian の最小値が既存手法を下回っており、既存手法のほうが高いメッシュ品質であることがわかります。以上の検討から、提案手法と既存手法は同程度の品質をもつメッシュ生成が可能であることがわかりました。
既存手法との比較:計算時間
三角形形状、円形状、五角形形状、リング形状の 4 つの領域に対し、四角形メッシュを抽出するために用いる三角形メッシュの節点数を増大させたときの計算時間の変化について評価します。計算環境は CPU Intel Core i7 (7567U)、メモリ (LPDDR3-2133) 16GB (Macbook pro) を利用しました。図 6 に、同一の計算を 5 回実施したうちの最も短い計算時間を示します。提案手法は既存手法に比べ、三角形メッシュの節点数が増大しても計算時間の増大は緩やかであることがわかりました。この結果により、提案手法はメッシュ抽出に必要な計算量を低減できることを確認しました。
図 6:提案手法と既存手法の計算時間の関係。四角形メッシュを抽出する元になる三角形メッシュの自由度が増大すると計算時間が増大する。提案手法は計算時間の増加が緩やかである。
その他手法との比較:メッシュ品質
三角形形状、円形状、五角形形状、リング形状の 4 つの領域に対し、以下の手法によって得られたメッシュの品質を比較しました。ここで、既存手法として gmsh の frontal delaunay quad meshing 法を利用しました。図 7 に、それぞれのメッシュ分割例を示します。
- 提案手法によって得られた四角形メッシュ
- Delaunay 法によって得られた四角形メッシュ
- 手作業で補助線を付与して得られた四角形メッシュ
表 3 に、メッシュ分割例の aspect ratio の最大値、平均値、最小値を示します。提案手法と既存手法を比べると、提案手法は aspect ratio の最大値を低減できており、品質の良いメッシュ分割が得られたことがわかります。提案手法と手作業をの例を比べると、提案手法は五角形領域以外の例で aspect ratio の最大値を低減できており、品質の良いメッシュ分割が得られたことがわかります。五角形領域では同程度の品質をもつメッシュ分割が得られました。
表 4 に、メッシュ分割例の scaled Jacobian の最大値、平均値、最小値を示します。提案手法と既存手法を比べると、提案手法は scaled Jacobian の最小値を増大できており、品質の良いメッシュ分割が得られたことがわかります。提案手法と手作業の例を比べると、提案手法は五角形領域以外の例で scaled Jacobian の最小値を増大できており、品質の良いメッシュ分割が得られたことがわかります。五角形領域では同程度の品質をもつメッシュ分割が得られました。
以上の数値例から基本形状の範囲内では、提案手法は既存手法や手作業によるメッシュ分割と同程度以上のメッシュ品質を有する完全四角形メッシュ分割が可能であることがわかりました。
おわりに
この記事では、有限要素法解析の根幹をなすメッシュ生成技術に関連して、frame field を用いた二次元メッシュ生成について紹介しました。計算精度や計算時間の観点から、品質のよい四角形メッシュの利用が促進できるよう、今後も継続した研究開発を予定しています。また、本研究開発の素直な発展である、三次元モデルにおける完全六面体メッシュ生成技術の開発にも取り組んでいます。また、本記事で紹介したメッシュ生成技術に関してのご相談も随時お受けしております。
また RICOS ではメッシュ生成技術の開発とあわせて、FEM と幾何学の機械学習に関する研究を進めています。有限要素メッシュを精度よく自動生成する技術について研究開発を進めることで、機械学習用途として、人手を介さず精度のよい膨大な解析データを自動的に生成することが可能となります。現在これら研究開発を促進するため、当該分野に関する研究する研究者、エンジニアを募集しております。この記事を読んで少しでも興味を持たれましたら、RICOS への応募を検討してみて下さい。
- Bommes, D., Zimmer, H. and Kobbelt, L., Mixed-integer quadrangulation, ACM SIGGRAPH 2009, ACM (2009), pp. 1–10.
- Bommes, D., Zimmer, H. and Kobbelt, L., Practical mixed-integer optimization for geometry processing, International Conference on Curves and Surfaces, Springer (2010), pp. 193–206.
- Bommes, D., Levy, B., Pietroni, N., Puppo, E., Silva, C., Tarini, M. and Zorin, D., State of the art in quad meshing, Eurographics 33rd Annual Conference of the European Association for Computer Graphics (2012).
- Bommes, D., Campen, M., Ebke, H. C., Alliez, P. and Kobbelt, L., Integer-grid maps for reliable quad meshing, ACM Transactions on Graphics, Vol. 32, No. 4 (2013a), p. 98.
- Bommes, D., Levy, B., Pietroni, N., Puppo, E., Silva, C., Tarini, M. and Zorin, D., Quad-mesh generation and processing: A survey, In Computer Graphics Forum, Vol. 32, No. 6 (2013b), pp. 51–76.
- Ray, N., Li, W. C., Levy, B., Sheffer, A. and Alliez, P., Periodic global parameterization, ACM Transactions on Graphics, Vol. 25, No. 4 (2006), pp. 1460–1485.
- Ray, N., Dmitry, S. and Levy, B., Practical 3D frame field generation, ACM Transactions on Graphics, Vol. 35, No. 6 (2016), p. 233.
- Vaxman, A., Campen, M., Diamanti, O., Panozzo, D., Bommes, D., Hildebrandt, K. and Ben-Chen, M., Directional field synthesis, design, and processing, Computer Graphics Forum, Vol. 35, No. 2 (2016), pp. 545–572.
- 森田 直樹, 堀江 正信, 須藤 海, 有限要素離散化したframe fieldによる二次元メッシュ生成, 日本機械学会論文集, ISSN 2187-9761, https://doi.org/10.1299/transjsme.19-00337, https://www.jstage.jst.go.jp/article/transjsme/advpub/0/advpub_19-00337/_article/-char/ja, 2020.