特許 7136660 テトラメッシュ生成装置、及びプログラム

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】

(24)【登録日】2022-09-05

(45)【発行日】2022-09-13

(54)【発明の名称】テトラメッシュ生成装置、及びプログラム

(51)【国際特許分類】

   G06F 30/23 20200101AFI20220906BHJP

【ＦＩ】

G06F30/23

【請求項の数】 6

(21)【出願番号】P 2018205606

(22)【出願日】2018-10-31

(65)【公開番号】P2020071692

(43)【公開日】2020-05-07

【審査請求日】2021-06-16

(73)【特許権者】

【識別番号】000002107

【氏名又は名称】住友重機械工業株式会社

(74)【代理人】

【識別番号】100105887

【弁理士】

【氏名又は名称】来山 幹雄

(72)【発明者】

【氏名】大西 良孝

【審査官】合田 幸裕

(56)【参考文献】

【文献】特開平０７－２３０４８７（ＪＰ，Ａ）

【文献】特開２００７－０９４６２８（ＪＰ，Ａ）

【文献】国際公開第２０１３／０５８２１５（ＷＯ，Ａ１）

【文献】特開２００８－２５７５９１（ＪＰ，Ａ）

【文献】特開平１１－３３９０７３（ＪＰ，Ａ）

(58)【調査した分野】(Int.Cl.，ＤＢ名)

Ｇ０６Ｆ ３０／２３

ＩＥＥＥ Ｘｐｌｏｒｅ

ＪＳＴＰｌｕｓ（ＪＤｒｅａｍＩＩＩ）

(57)【特許請求の範囲】

【請求項1】

三次元形状を定義する形状定義データが入力される入力部と、
前記入力部に入力された前記形状定義データで定義される前記三次元形状のエッジ上、及び表面上に複数の粒子を配置し、前記三次元形状の表面から離れた内奥部に、内部分布規則に則って粒子を配置し、前記三次元形状の内奥部と表面との間の空間に、前記内部分布規則に則った位置からずれた位置に粒子を配置し、配置された粒子の中心を頂点とするテトラメッシュを生成する処理部と
を有するテトラメッシュ生成装置。

【請求項2】

前記処理部は、
前記三次元形状のエッジ上に粒子を配置した後、エッジ上の粒子を発生源として、前記三次元形状の表面上に、表面分布規則に則って粒子を順次配置し、
前記三次元形状の表面上の粒子を変位させて、異なる発生源から順次配置された粒子同士が衝突する境界における前記表面分布規則の破綻を緩和させる請求項１に記載のテトラメッシュ生成装置。

【請求項3】

前記処理部は、
前記三次元形状の内部に、前記内部分布規則に則って粒子を配置し、
前記三次元形状の内部の粒子を変位させて、前記表面分布規則の破綻を緩和させた後の表面上の粒子の集合と、前記三次元形状の内部の粒子の集合との界面に生じている粒子の分布密度の不連続を緩和させる請求項２に記載のテトラメッシュ生成装置。

【請求項4】

さらに、出力部を有し、
前記処理部は、
テトラメッシュを生成する処理において、前記三次元形状のエッジ上の粒子、前記表面分布規則の破綻を緩和させた後の表面上の粒子、及び粒子の分布密度の不連続を緩和させた後の前記三次元形状の内部の粒子の位置に基づいて、粒子の中心を頂点とするテトラメッシュを生成し、
生成したテトラメッシュを定義するテトラメッシュ定義データを前記出力部に出力する請求項３に記載のテトラメッシュ生成装置。

【請求項5】

前記処理部は、
粒子の中心を頂点とするテトラメッシュを生成した後、前記テトラメッシュ定義データを出力する前に、前記三次元形状の外側の不要な四面体を削除する請求項４に記載のテトラメッシュ生成装置。

【請求項6】

三次元形状を定義する形状定義データを入力させる機能と、
入力された前記形状定義データで定義される前記三次元形状のエッジ上、及び表面上に複数の粒子を配置し、前記三次元形状の表面から離れた内奥部に、内部分布規則に則って粒子を配置し、前記三次元形状の内奥部と表面との間の空間に、前記内部分布規則に則った位置からずれた位置に粒子を配置し、粒子の中心を頂点とするテトラメッシュを生成す
る機能と、
生成されたテトラメッシュを定義するテトラメッシュ定義データを出力する機能と
をコンピュータに実行させるプログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、テトラメッシュ生成装置、及びプログラムに関する。

【背景技術】

【0002】

テトラメッシュの生成方法が多数提案されている（特許文献１、非特許文献１参照）。これらのテトラメッシュの生成方法のほとんどは、有限要素法への適用を想定している。そのため、結晶格子構造を基にしている分子動力学法やくりこみ群分子動力学法（以下、両者をまとめて分子動力学法という場合がある。）とは相性が良くない場合がある。

【先行技術文献】

【特許文献】

【0003】

【文献】特開２００３－７５５２１号公報

【非特許文献】

【0004】

【文献】Parikh, P., Pizadeh, S., and Lohner, R., “Package for 3-D Unstructured Grid Generation, Finite-Element Flow Solution and Flow Field Visualization”, NASA CR-182090, 1990

【発明の概要】

【発明が解決しようとする課題】

【0005】

従来の方法で生成したテトラメッシュを用いて分子動力学法による解析を行うと、解析精度が著しく低下する場合がある。特に、直方体を分割して作られたテトラメッシュ（特許文献１参照）はせん断方向への剛性が低いため、このメッシュを用いて分子動力学法による解析を行うと、解析が破綻する場合がある。

【0006】

本発明の目的は、分子動力学法やくりこみ群分子動力学法による解析に適したテトラメッシュ生成装置、及びプログラムを提供することである。

【課題を解決するための手段】

【0007】

本発明の一観点によると、
三次元形状を定義する形状定義データが入力される入力部と、
前記入力部に入力された前記形状定義データで定義される前記三次元形状のエッジ上、及び表面上に複数の粒子を配置し、前記三次元形状の表面から離れた内奥部に、内部分布規則に則って粒子を配置し、前記三次元形状の内奥部と表面との間の空間に、前記内部分布規則に則った位置からずれた位置に粒子を配置し、配置された粒子の中心を頂点とするテトラメッシュを生成する処理部と
を有するテトラメッシュ生成装置が提供される。

【0009】

本発明の他の観点によると、
三次元形状を定義する形状定義データを入力させる機能と、
入力された前記形状定義データで定義される前記三次元形状のエッジ上、及び表面上に複数の粒子を配置し、前記三次元形状の表面から離れた内奥部に、内部分布規則に則って粒子を配置し、前記三次元形状の内奥部と表面との間の空間に、前記内部分布規則に則った位置からずれた位置に粒子を配置し、粒子の中心を頂点とするテトラメッシュを生成する機能と、
生成されたテトラメッシュを定義するテトラメッシュ定義データを出力する機能と
をコンピュータに実行させるプログラムが提供される。

【発明の効果】

【0010】

内奥部に、シミュレーション対象の原子配置を反映した形状のテトラメッシュを生成することが可能である。内奥部と表面との間の空間においては、内部分布規則に則った位置からずれて粒子を配置するため、内奥部のテトラメッシュと、表面のテトラメッシュとのメッシュ形状の不連続性を緩和させることができる。

【0011】

このテトラメッシュの節点に粒子を配置して、分子動力学法またはくりこみ群分子動力学法により粒子の挙動をシミュレーションすることにより、シミュレーションの精度を高めることができる。

【図面の簡単な説明】

【0012】

【図1】図１は、実施例によるテトラメッシュ生成装置のブロック図である。

【図2】図２は、実施例によるテトラメッシュ生成方法のフローチャートである。

【図3】図３は、ＳＴＬ形式で定義された三次元形状の一例を示す斜視図である。

【図4】図４は、複数の区画を階調の異なるグレーレベルで色分けした三次元形状の斜視図である。

【図5】図５Ａは、区画ごとに階調の異なるグレーレベルで色分けした三次元形状の一部分を拡大した図であり、図５Ｂは、エッジ上の三角形要素の節点に球形の粒子を表示した図である。

【図6】図６は、エッジ上に配置した複数の節点を示す斜視図である。

【図7】図７は、三次元形状の表面上に粒子を配置する表面分布規則を説明するための図である。

【図8】図８Ａは、三次元形状の表面を覆いつくした粒子の分布の一例を示す図であり、図８Ｂは、図８Ａの一部分を拡大した図である。

【図9】図９は、三次元形状の表面から逸脱した粒子を表面上に戻す処理を説明するための模式図である。

【図10】図１０は、三次元形状が配置された空間を複数のセルに分割した状態の模式図である。

【図11】図１１Ａ～図１１Ｃは、表面分布規則の破綻を緩和する処理（ステップＳ５）を行う前の粒子の分布を示す模式図であり、図１１Ｄ～図１１Ｆは、それぞれ図１１Ａ～図１１Ｃに示した例の緩和処理後の粒子の分布を示す模式図である。

【図12】図１２は、三次元形状の表面を構成しない余分な三角形の一例を示す斜視図である。

【図13】図１３は、ステップＳ５の緩和処理を行うときの運動方程式のタイムステップ数を変えて得られた表面メッシュの歪み量Ｑを示す箱ひげ図である。

【図14】図１４Ａは、三次元形状の内部に配置した粒子、表面上の粒子、及びエッジ上の粒子の分布の一例を示す模式図であり、図１４Ｂは、粒子の分布密度の不連続を緩和させる処理を行った後の粒子の分布を示す模式図である。

【図15】図１５Ａは、表面メッシュを構成する三角形がテトラメッシュの四面体の面として採用されていない例を示す模式図であり、図１５Ｂは、四面体の頂点の組み合わせを調整した後のテトラメッシュを示す模式図である。

【図16】図１６Ａは、最終的に作成されたテトラメッシュの斜視図であり、図１６Ｂは、余分な四面体を削除する前のテトラメッシュの斜視図である。

【図17】図１７は、三次元形状の内部に形成されたテトラメッシュの形状を示す断面図、及び一部の断面の拡大図である。

【発明を実施するための形態】

【0013】

図１～図１６を参照しながら、実施例によるテトラメッシュ生成装置及びテトラメッシュ生成方法について説明する。

【0014】

図１は、実施例によるテトラメッシュ生成装置のブロック図である。実施例によるテトラメッシュ生成装置は、入力部２０、処理部２１、出力部２２、及び記憶部２３を含む。入力部２０から、テトラメッシュ生成の条件、例えば三次元形状を定義する形状定義データ、メッシュの節点間隔、節点の配置の基本的な条件等が入力される。さらに、オペレータから入力部２０に各種指令（コマンド）等が入力される。入力部２０は、例えば通信装置、リムーバブルメディア読取装置、キーボード等で構成される。

【0015】

処理部２１は、入力された条件に基づいてテトラメッシュを生成する処理を行い、処理結果を出力部２２に出力する。処理結果には、テトラメッシュの形状を定義するテトラメッシュ定義データが含まれる。処理部２１は、例えばコンピュータを含み、テトラメッシュを生成する種々の機能をコンピュータに実行させるためのプログラムが記憶部２３に記憶されている。出力部２２は、通信装置、リムーバブルメディア書込み装置、ディスプレイ等を含む。

【0016】

図２は、実施例によるテトラメッシュ生成方法のフローチャートである。
ステップＳ１において、処理部２１が、入力部２０から入力された三次元形状を定義する形状定義データを取得する。形状定義データは、例えばＳＴＬ（Standard Triangulated Language）形式のデータであり、三次元ＣＡＤを用いて作成される。ＳＴＬ形式のデータは、小さな三角形要素の集合体で三次元形状の表面を再現する。

【0017】

図３は、ＳＴＬ形式で定義された三次元形状３０の一例を示す斜視図である。図３に示した三次元形状のモデルは混錬機のローラーである。図３では、三次元形状３０の表面に配置される複数の三角形要素の頂点に相当する位置に球形の粒子３１を示している。平坦な領域及び曲率半径の大きな領域には、大きな三角形要素が配置され、曲率半径の小さな領域には、小さな三角形要素が配置される。

【0018】

ステップＳ２（図２）において、処理部２１は、形状定義データに基づいて、三次元形状３０のエッジを検出する。以下、エッジを検出する方法の一例について説明する。

【0019】

形状定義データの三角形要素の各々に、法線ベクトルを定義することができる。隣り合う三角形要素の法線ベクトルのなす角度が所定のしきい値以下であれば、この２つの三角形要素を同一の区画に属させる。隣り合う三角形要素の法線ベクトルのなす角度が所定のしきい値より大きい場合、この２つの三角形要素を異なる区画に属させる。このしきい値は、例えば３０度程度とする。区画の境界をエッジと定義する。

【0020】

図４は、複数の区画を階調の異なるグレーレベルで色分けした三次元形状３０の斜視図である。図５Ａは、複数の区画を階調の異なるグレーレベルで色分けした三次元形状３０の一部分を拡大した図である。図４及び図５Ａに示したように、区画の境界にエッジ３２が定義される。図５Ｂは、エッジ３２上の三角形要素の節点に球形の粒子３１を表示した図である。

【0021】

ステップＳ３（図２）において、エッジ３２上に等間隔で複数の節点を配置する。この節点は、テトラメッシュの節点に相当する。節点の間隔は、テトラメッシュを生成する条件として、入力部２０（図１）から入力される。

【0022】

図６は、エッジ３２上に配置した複数の節点を示す斜視図である。図６においては、節点に相当する位置に、球形の粒子３５を表示している。本実施例で生成したテトラメッシュを分子動力学法に適用する場合には、テトラメッシュの節点の位置に粒子を配置するため、本明細書において、テトラメッシュの節点を「粒子」という場合がある。テトラメッシュの節点の位置は、粒子の中心位置に相当する。

【0023】

ステップＳ４（図２）において、エッジ３２上の粒子３５（図５）を発生源として、三次元形状３０（図３）の表面に、表面分布規則に則って粒子３６を順次配置する。図７を参照して、表面分布規則について説明する。

【0024】

図７は、三次元形状３０（図３）の表面上に粒子３６を配置する表面分布規則を説明するための図である。ステップＳ３（図２）において、エッジ３２上に複数の粒子３５が等間隔で配置されている。エッジ３２上の複数の粒子３５のうち隣り合う２つの粒子３５ａと３５ｂとを発生源として、２つの発生源から等距離（エッジ３２上の粒子３５の間隔）で、かつ三次元形状３０の三角形要素上の位置に、新たに粒子３６ａを配置する。続いて、粒子３５ａと３６ａとに基づいて、同様に新たな粒子３６ｂを配置する。

【0025】

このように、エッジ３２上の粒子３５、または既に配置した表面上の粒子３６の位置に基づいて、三次元形状３０の表面を覆いつくすまで同様の処理を繰り返す。この処理により、発生源となるエッジ３２上の粒子３５から、表面上の粒子３６が配置された領域が広がる（成長する）。新たに配置される粒子３７の近傍に、既に粒子３５または３６が配置されている場合には、その粒子３７を配置しない。ここで、「近傍」とは、エッジ３２上の粒子３５の間隔よりも近い範囲を意味する。共通の発生源から成長した領域を均質領域ということとする。

【0026】

図８Ａは、三次元形状３０の表面を覆いつくした粒子３６の分布の一例を示す図である。図８Ｂは、図８Ａの一部分を拡大した図である。異なる発生源から成長した２つの均質領域が衝突する位置に欠陥や粒界等が発生した境界３８が現れる。この境界３８においては、表面分布規則が破綻している。

【0027】

ステップＳ５（図２）において、表面上の粒子３６を変位させて、均質領域の境界における表面分布規則の破綻を緩和させる。ここで、破綻を緩和させるとは、境界における粒子配置の規則性の乱れを緩和させることを意味する。以下、表面分布規則の破綻を緩和させる処理の一例について説明する。

【0028】

エッジ３２上の粒子３５及び表面上の粒子３６の間の相互作用ポテンシャルを定義する。この相互作用ポテンシャルを用いて運動方程式を解くことにより、仮想的な粒子を移動させる。なお、エッジ３２上の粒子３５は移動しないものとする。

【0029】

相互作用ポテンシャルとしてバネモデルを用い、斥力のみを考慮する。また、エッジ３２上の粒子以外の粒子の粒径をゼロとし、時間の経過とともに粒子を膨張させる、即ち粒径を大きくする。粒子に作用する圧力が所定のしきい値を超えたら、粒子の膨張を停止させる。例えば、隣の粒子と接触するまで膨張したら、当該粒子の膨張を停止させる。粒子同士がほぼ接触した状態で運動方程式を解く場合には、運動方程式のタイムステップを繰り返しても粒子がほとんど移動しない。粒子の粒径の初期値をゼロに設定することにより、粒子を容易に変位させることができる。

【0030】

粒子ｉに作用する圧力Ｐｉは、以下の式で定義される。

【数1】

ここで、ｊは、粒子ｉと相互作用する粒子を意味し、Ｖ_ｉは粒子ｉの体積であり、ｒ_ｉｊは粒子ｉと粒子ｊとの距離であり、φ（ｒ_ｉｊ）は粒子ｉと粒子ｊとの間のポテンシャル関数である。

【0031】

運動方程式を解いて三次元形状３０（図３）の表面上の粒子を移動させると、粒子が表面から逸脱してしまう場合がある。表面から逸脱した粒子は、三次元形状３０の表面上に戻す。次に、粒子を表面上に戻す処理について、図９を参照して説明する。

【0032】

図９は、三次元形状の表面から逸脱した粒子を表面上に戻す処理を説明するための模式図である。表面３９上の粒子４０が、運動方程式を解いて１タイムステップ進めると粒子４０ａの位置まで移動する。粒子４０ａは三次元形状の表面３９から逸脱している。このとき、粒子４０ａを三次元形状の表面３９の三角形要素が作る平面に垂直投影した投影像が三角形要素の内部にあり、かつ粒子４０ａに最も近い投影像４０ｂを定義する。この投影像４０ｂの位置を、運動方程式を解く前の粒子４０が移動した後の位置と定義する。通常は、移動後の粒子４０ａが表面３９から逸脱する距離は僅かであるため、粒子４０ａの位置に最も近い三角形要素が作る平面上に垂直投影した投影像が、上述の投影像４０ｂとして採用される。

【0033】

次に、図１０を参照して、表面から逸脱した粒子４０ａ（図９）の位置に最も近い三角形要素を見つる方法の一例について説明する。

【0034】

図１０は、三次元形状３０（図３）が配置された空間を複数のセル４３に分割した状態の模式図である。図１０では、三次元形状３０を、その表面に対して垂直方向から見た模式図である。三次元形状３０が配置された空間が複数のセル４３に分割されている。予め、三角形要素と、その三角形要素が交差するセル４３とを関連付けて登録しておく。さらに、いずれの三角形要素とも交差していないセル４３については、当該セル４３に最も近い位置にある三角形要素と関連付けておく。セル４３と三角形要素とを関連付けたデータは、記憶部２３（図１）に記憶されている。図１０に示した例では、ハッチングが付されたセル４３が、実線で示した三角形要素４４に関連付けられている。

【0035】

運動方程式を解いた後、表面から逸脱した粒子４０ａ（図９）が含まれるセル４３を求める。セル４３と三角形要素とを関連付けたデータに基づいて、このセル４３に関連付けられている三角形要素を抽出する。このように、三次元形状３０が配置されている空間を複数のセル４３に分割して、セル４３と三角形要素と関連付けておくことにより、表面から逸脱した粒子４０ａ（図９）の位置に最も近い三角形要素を見つけることができる。

【0036】

図１１Ａ～図１１Ｃは、表面分布規則の破綻を緩和する処理（ステップＳ５）を行う前の粒子の分布を示す模式図である。図１１Ａに示した例では、粒子の分布密度が相対的に高い領域と、相対的に低い領域とが観察される。図１１Ｂに示した例では、粒子の分布密度が周囲よりも高くなっている領域が観察される。図１１Ｃに示した例では、粒子の分布密度が局所的に高くなっている領域が観察される。

【0037】

図１１Ｄ～図１１Ｆは、それぞれ図１１Ａ～図１１Ｃに示した例の緩和処理後の粒子の分布を示す模式図である。いずれの例においても、粒子の分布が均等な分布に近づいていることがわかる。このように、ステップＳ５の緩和処理を行うことにより、三次元形状３０の表面上に配置された粒子の分布を均等な分布に近付けることができる。

【0038】

ステップＳ６（図２）において、緩和処理を行った後の粒子の位置に基づいて、表面メッシュを作成する。以下、表面メッシュの作成方法について説明する。

【0039】

緩和処理を行った後の粒子の位置に基づいて、ドロネー三角形分割を行う。具体的には、まず、３個の粒子を選択し、この３個の粒子の中心を表面に含む球を考える。ただし、３個の粒子の中心位置から得られる平面と球との交円は大円であるとする。球の内部に他の粒子が存在しなければ、選択した３個の粒子の中心位置を頂点とする三角形を、表面メッシュの一部として採用する。

【0040】

ただし、この方法では、実際には三次元形状３０の表面を構成していない三角形が表面メッシュの一部として採用されてしまう場合がある。このような三角形を削除することにより、表面メッシュを完成する。

【0041】

図１２は、三次元形状３０の表面を構成しない余分な三角形の一例を示す斜視図である。三次元形状３０の表面４５と表面４６とが交線４７で交わっている。このように２つの表面が交わっている交線４７上の粒子を含む余分な三角形４８が表面メッシュの一部として採用されやすい。三次元形状３０の表面メッシュを構成する三角形の各々は、３つの辺を介して３つの三角形に連続している。ところが、余分な三角形４８は、２つの辺を介してそれぞれ隣の三角形に連続しているが、残りの１つの辺を介して他の三角形に連続していない。この辺は、他の三角形と共有されていない。このように、辺を介して隣り合う三角形の個数が３個ではない三角形を余分な三角形として削除する。

【0042】

このようにして作成した表面メッシュの歪み量Ｑを計算により求めた。歪み量Ｑは、以下の式で定義される。

【数2】

ここで、θ_ｍａｘ及びθ_ｍｉｎは、それぞれ表面メッシュを構成する三角形の内角の最大値及び最小値である。θ_ｅは正三角形の１つの内角の角度、すなわち６０°である。表面メッシュのすべての三角形が正三角形であるとき歪み量Ｑ＝０になり、三角形のすべての頂点が１本の直線上に並ぶと歪み量Ｑ＝１になる。

【0043】

図１３は、ステップＳ５の緩和処理を行うときの運動方程式のタイムステップ数を変えて得られた緩和処理後の表面メッシュの歪み量Ｑを示す箱ひげ図である。参考のために、市販のメッシュ生成ソフトを用いて生成した表面メッシュの歪み量Ｑを比較例として示す。

【0044】

緩和処理のタイムステップ数を増やすと、歪み量Ｑが低下することがわかる。また、本実施例による方法で作成した表面メッシュの歪み量Ｑは、従来の市販のソフトを使用して作成した表面メッシュの歪み量Ｑより小さいことがわかる。テトラメッシュを用いて分子動力学法により解析を行う場合には、経験的に歪み量Ｑが０．５以下であれば十分であることがわかっている。本実施例による方法により、この条件を満たす表面メッシュを生成することができる。

【0045】

ステップＳ７（図２）において、三次元形状３０の内部に、内部分布規則に則って粒子を配置する。内部分布規則として、例えば、面心立方格子の原子位置に粒子を配置するという規則を採用する。この規則は、面心立方格子構造の金属を分子動力学法によってシミュレーションする場合に特に有効である。なお、体心立方格子構造や六方最密充填構造の金属のシミュレーションを行う場合には、内部分布規則としてこれらの構造の原子位置に粒子を配置するという規則を採用するとよい。

【0046】

次に、粒子の配置方法についてより具体的に説明する。まず、三次元形状３０（図３）の内部、及び表面メッシュをすべて覆うように粒子を配置する。粒子の半径は、内部分布規則に適用される格子構造の原子半径とする。例えば、内部分布規則として面心立方格子を採用する場合、粒子の半径ｒは格子定数ａを用いて以下の式で定義される。

【数3】

【0047】

次に、表面上及びエッジ上の粒子３５、３６（図６、図７）と重なる粒子を削除する。なお、表面上及びエッジ上の粒子３５、３６の半径は、内部分布規則に則って配置される粒子の半径と同一である。内部分布規則に則って配置した粒子の中心間の距離から、結合している粒子対を求める。結合している粒子同士を同一のグループに分類する。表面上及びエッジ上の粒子３５、３６と重なる粒子を削除しているため、内部分布規則に則って配置された粒子は２つのグループに分類される。２つのグループの粒子のうち、三次元形状３０の外側に配置されていると判定されたグループの粒子を削除する。この結果、三次元形状３０の内部のみに粒子を配置することができる。

【0048】

図１４Ａは、三次元形状３０の内部に配置した粒子５２、表面上の粒子３６、及びエッジ上の粒子３５の分布の一例を示す模式図である。内部分布規則に則って配置した粒子のうち、表面上の粒子３６及びエッジ上の粒子３５と重なる粒子を削除しているため、内部の粒子５２の集合と表面上及びエッジ上の粒子３６、３５の集合との界面に、粒子５２の直径より小さい寸法の初期隙間５３が生じる。すなわち、内部の粒子５２の集合と、表面上及びエッジ上の粒子３６、３５の集合との界面に、粒子の分布密度が不連続に変化する部分が発生する。

【0049】

ステップＳ８（図２）において、内部の粒子５２を変位させて粒子の分布密度の不連続を緩和させる。内部の粒子５２を変位させる方法は、ステップＳ５で表面上の粒子を変位させた方法と同一である。なお、内部の粒子５２を変位させる処理では、表面から逸脱した粒子を表面上に戻す処理は不要である。また、表面上及びエッジ上の粒子３６、３５は移動させない。

【0050】

図１４Ｂは、粒子の分布密度の不連続を緩和させる処理を行った後の粒子の分布を示す模式図である。内部の粒子５２の集合と表面上及びエッジ上の粒子３６、３５の集合との界面に生じていた初期隙間５３（図１４Ａ）が、ほぼ消滅していることがわかる。

【0051】

ステップＳ９（図２）において、内部の粒子５２、表面上の粒子３６、及びエッジ上の粒子３５（図１４Ｂ）の中心を頂点とするテトラメッシュを作成する。テトラメッシュの作成には、例えばドロネー三角形分割手法を用いることができる。ただし、この手法により作成したテトラメッシュは、表面メッシュを構成する三角形がテトラメッシュの四面体の面として採用されない状況が生じ得る。

【0052】

図１５Ａは、表面メッシュを構成する三角形がテトラメッシュの四面体の面として採用されていない例を示す模式図である。三次元形状３０（図３）の表面上に節点Ｉ、Ｊ、Ｋが配置され、三次元形状３０の内部に節点Ｌが配置され、節点Ｉ、Ｊ、Ｋを含む表面の外側に節点Ｍが配置されている。節点Ｉ、Ｋ、Ｌ、Ｍを頂点とする１つの四面体ＩＫＬＭが作成され、節点Ｊ、Ｋ、Ｌ、Ｍを頂点とする他の１つの四面体ＪＫＬＭが作成されている。図１５Ａにおいて、２つの四面体の界面ＫＬＭにハッチングを付している。この場合、節点Ｉ、Ｊ、Ｋを頂点とする三角形ＩＪＫが四面体の１つの面として採用されていない。四面体ＩＫＬＭ及び四面体ＪＫＬＭは、表面メッシュの三角形ＩＪＫを貫いている。

【0053】

ステップＳ１０（図２）において、表面メッシュを構成するすべての三角形がテトラメッシュの四面体の面を構成するように、四面体の頂点の組み合わせを調整する。

【0054】

図１５Ｂは、四面体の頂点の組み合わせを調整した後のテトラメッシュを示す模式図である。以下、四面体の頂点の組み合わせの調整方法について説明する。まず、表面メッシュを貫いている四面体を抽出する。図１５Ａの例では、四面体ＩＫＬＭ及び四面体ＪＫＬＭが抽出される。この２つの四面体の頂点の組み合わせを変えて、２つの四面体ＩＪＫＬ及び四面体ＩＪＫＭを作成する。図１５Ｂにおいて、２つの四面体の界面ＩＪＫにハッチングを付している。図１５Ｂに示した例では、表面メッシを構成する三角形ＩＪＫが、テトラメッシュの四面体の面として採用されている。

【0055】

３つの四面体が、表面メッシュの隣り合う２つの三角形を貫いている場合も、同様に３つの四面体の頂点の組み合わせを調整する。表面メッシュの隣り合う３つ以上の三角形を複数の四面体が貫いている場合も、同様の処理を行う。

【0056】

ステップＳ１１（図２）において、不要な四面体、すなわち、三次元形状３０の外側の四面体を削除する。次に、不要な四面体を削除する方法について説明する。

【0057】

まず、作成されたテトラメッシュの四面体から表面に位置する四面体を抽出する。例えば、四面体の１つの面を共有する他の四面体が存在しない場合、その面は表面であると判定する。表面と判定された面が表面メッシュの三角形と一致しない場合、当該四面体を削除する。四面体を削除すると、新たに表面が現れるため、新たな表面と表面メッシュの三角形とを比較する。四面体を削除する処理を繰り返すことにより、三次元形状３０の外側の四面体が削除され、内部の四面体のみが残る。

【0058】

図１６Ａは、最終的に作成されたテトラメッシュの斜視図である。図１６Ｂは、余分な四面体を削除する前のテトラメッシュの斜視図である。図１６Ｂに示したテトラメッシュは、表面の窪みを埋め込むように三次元形状３０の外側に配置された四面体を含んでいる。三次元形状３０の外側の四面体を削除することにより、図１６Ａに示したテトラメッシュが得られる。

【0059】

図１７は、三次元形状３０の内部に形成されたテトラメッシュの形状を示す破断斜視図、及び断面の一部を拡大した断面図である。三次元形状３０の表面から離れた内奥部においては、ほぼ面心立方格子に相当するテトラメッシュが形成されていることがわかる。ステップＳ８（図２）での緩和処理では、内奥部の粒子は実質的にほとんど移動しない。このため、内奥部には、内部分布規則に則った粒子５２の配置が維持されている。内奥部と表面との間の空間、特に表層部においては、粒子５２が変位することにより、内部分布規則に則った位置からややずれた位置に粒子５２が配置されることになる。

【0060】

次に、上記実施例の優れた効果について説明する。
上記実施例では、三次元形状の内奥部に、シミュレーション対象の原子配置を反映した形状のテトラメッシュを生成することができる。また、三次元形状の表面やエッジには、ほぼ均等に粒子を配置することができる。内奥部と表面との間の空間においては、内部分布規則に則った位置からずれて粒子を配置するため、内奥部のテトラメッシュと、表面のテトラメッシュとのメッシュ形状の不連続性を緩和させることができる。

【0061】

このテトラメッシュの節点に粒子を配置して、分子動力学法またはくりこみ群分子動力学法を用いて粒子の挙動をシミュレーションすることにより、シミュレーションの精度を高めることができる。

【0062】

本発明は上述の実施例に制限されるものではない。例えば、種々の変更、改良、組み合わせ等が可能なことは当業者に自明であろう。

【符号の説明】

【0063】

２０ 入力部
２１ 処理部
２２ 出力部
２３ 記憶部
３０ 三次元形状
３１ 三次元形状を反映する三角形要素の頂点に位置する粒子
３２ 三次元形状のエッジ
３５、３５ａ、３５ｂ エッジ上の粒子
３６、３６ａ、３６ｂ 表面上の粒子
３７ 配置しない粒子
３８ 境界
３９ 三次元形状の表面
４０ 運動方程式を解く前の粒子
４０ａ 運動方程式を解いて１タイムステップ進ませた後の粒子
４０ｂ 表面への垂直投影像
４３ セル
４４ 三次元形状の表面の三角形要素
４５、４６ 三次元形状の表面
４７ 交線
４８ 表面を構成しない余分な三角形
５２ 内部の粒子
５３ 初期隙間

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】

【図12】

【図13】

【図14】

【図15】

【図16】

【図17】