(57)【要約】  第１の物体を表現した複数の粒子を含む第１の粒子群の物理量に関するデータに基づき、前記第１の粒子群に含まれる前記複数の粒子の間を接続する接続要素を生成する生成部と、第２の物体に含まれる第２の粒子の物理量に関するデータを参照し、前記第２の粒子と前記生成された接続要素との最短距離を算出して、前記第１の物体と前記第２の粒子との間に作用する反発力を算出する算出部と、前記算出された反発力を用いて、前記第１の物体及び前記第２の物体の運動を解析する解析部と、を有することを特徴とする物体運動解析装置。
2013042234.tif 000002
【特許請求の範囲】
【請求項１】
第１の物体を表現した複数の粒子を含む第１の粒子群の物理量に関するデータに基づき、前記第１の粒子群に含まれる前記複数の粒子の間を接続する接続要素を生成する生成部と、
第２の物体に含まれる第２の粒子の物理量に関するデータを参照し、前記第２の粒子と前記生成された接続要素との最短距離を算出して、前記第１の物体と前記第２の粒子との間に作用する反発力を算出する算出部と、
前記算出された反発力を用いて、前記第１の物体及び前記第２の物体の運動を解析する解析部と、
を有することを特徴とする物体運動解析装置。
【請求項２】
前記物体運動解析装置において、
前記算出部は、前記粒子と前記接続要素との最短距離が減少してゼロに近づくにつれて前記反発力が大きくなり、前記第２の粒子と前記接続要素との最短距離が増加して所定距離に近づくにつれて前記反発力がゼロに近づくように、前記反発力を算出することを特徴とする、
請求項１記載の物体運動解析装置。
【請求項３】
前記物体運動解析装置において、
構造物を表現するデータを分割することにより、前記複数の粒子に関するデータを生成する分割部を備えることを特徴とする、
請求項１記載の物体運動解析装置。
【請求項４】
前記物体運動解析装置において、
前記第１の物体及び前記第２の物体は、同一平面上に存在する物体であり、
前記接続要素は、前記複数の第１の粒子間を接続する線分要素であることを特徴とする、
請求項１記載の物体運動解析装置。
【請求項５】
前記物体運動解析装置において、
前記第１の物体及び前記第２の物体は、三次元空間に存在する物体であり、
前記接続要素は、前記複数の第１の粒子から選択された三個の第１の粒子を含む閉領域であることを特徴とする、
請求項１記載の物体運動解析装置。
【請求項６】
物体の運動を解析する物体運動解析方法において、
コンピュータが、
第１の物体に含まれる複数の第１の粒子の物理量に関するデータに基づき、前記複数の第１の粒子間を接続する接続要素を生成し、
第２の物体に含まれる第２の粒子の物理量に関するデータを参照し、前記第２の粒子と前記生成された接続要素との最短距離を算出して、前記第１の物体と前記第２の粒子との間に作用する反発力を算出し、
前記算出された反発力を用いて、前記第１の物体及び前記第２の物体の運動を解析する、
ことを特徴とする物体運動解析方法。
【請求項７】
前記物体運動解析方法において、
前記最短距離が減少してゼロに近づくと前記反発力が大きくなり、前記第２の粒子と前記接続要素との最短距離が増加して所定距離に近づくと前記反発力がゼロに近づくように、前記反発力を算出することを特徴とする、
請求項６記載の物体運動解析方法。
【請求項８】
前記物体運動解析方法はさらに、
前記コンピュータが、
構造物を表現するデータを分割することにより、前記複数の粒子に関するデータを生成することを特徴とする、
請求項６記載の物体運動解析方法。
【請求項９】
前記物体運動解析方法において、
前記第１の物体及び前記第２の物体は、同一平面上に存在する物体であり、
前記接続要素は、前記複数の第１の粒子間を接続する線分要素であることを特徴とする、
請求項６記載の物体運動解析方法。
【請求項１０】
前記物体運動解析方法において、
前記第１の物体及び前記第２の物体は、三次元空間に存在する物体であり、
前記接続要素は、前記複数の第１の粒子から選択された三個の第１の粒子を含む閉領域であることを特徴とする、
請求項６記載の物体運動解析方法。
【請求項１１】
物体の運動を解析する物体運動解析プログラムにおいて、
コンピュータに、
第１の物体に含まれる複数の第１の粒子の物理量に関するデータに基づき、前記複数の第１の粒子間を接続する接続要素を生成させ、
第２の物体に含まれる第２の粒子の物理量に関するデータを参照し、前記第２の粒子と前記生成された接続要素との最短距離を算出して、前記第１の物体と前記第２の粒子との間に作用する反発力を算出させ、
前記算出された反発力を用いて、前記第１の物体及び前記第２の物体の運動を解析させる、
ことを特徴とする物体運動解析プログラム。
【請求項１２】
前記物体運動解析プログラムにおいて、
前記最短距離が減少してゼロに近づくと前記反発力が大きくなり、前記第２の粒子と前記接続要素との最短距離が増加して所定距離に近づくと前記反発力がゼロに近づくように、前記反発力を算出させることを特徴とする、
請求項１１記載の物体運動解析プログラム。
【請求項１３】
前記物体運動解析プログラムはさらに、
前記コンピュータに、
構造物を表現するデータを分割することにより、前記複数の粒子に関するデータを生成させることを特徴とする、
請求項１１記載の物体運動解析プログラム。
【請求項１４】
前記物体運動解析プログラムにおいて、
前記第１の物体及び前記第２の物体は、同一平面上に存在する物体であり、
前記接続要素は、前記複数の第１の粒子間を接続する線分要素であることを特徴とする、
請求項１１記載の物体運動解析プログラム。
【請求項１５】
前記物体運動解析プログラムにおいて、
前記第１の物体及び前記第２の物体は、三次元空間に存在する物体であり、
前記接続要素は、前記複数の第１の粒子から選択された三個の第１の粒子を含む閉領域であることを特徴とする、
請求項１１記載の物体運動解析プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、物体同士の反発力を算出して物体の運動を解析する物体運動解析装置、物体運動解析方法及び物体解析運動プログラムに関する。
【背景技術】
【０００２】
従来、コンピュータ上で物体の運動を再現して挙動を解析する試みが、種々の分野で行われている。コンピュータによる物体の運動の再現結果は、例えば、クラブで打たれたゴルフボールの動き、津波を浴びた建築物の耐性、液体等を収容する容器に取り付けられたゴム膜の挙動等を調べるのに、広く利用することができる。【０００３】
このような運動解析対象となる物体には、流体や弾性体が含まれる。従来、流体や弾性体等の連続体を解く数値計算手法として、格子をベースにして微分方程式の近似解を求解する有限差分法や有限要素法、有限体積法等が多く用いられてきた。近年では、数値計算をＣＡＥ（Computer  Aided  Engineering）などの応用分野で活用するため、これらの数値計算手法も発展し、流体と構造物が相互作用する問題が解かれるようになってきた。しかしながら、格子を用いる手法においては、自由表面などの界面の存在する問題や、流体-構造連成問題等の移動境界が発生する場合には、取り扱いが複雑なため、プログラム作成が困難となる場合が多い。【０００４】
これに対して、物体を粒子の集合として解析する粒子法が知られている。代表的な粒子法としては、ＭＰＳ（Moving  Particle  Semi-implicit）法やＳＰＨ（Smoothed  Particle  Hydrodynamics）法等が知られている（例えば、非特許文献１、２参照）。粒子法では、移動境界の取り扱いに特別な処置を必要としないため、粒子法は流体−構造連成問題にも広く用いられるようになっている。【０００５】
これに関連し、流体を表した第１の粒子毎に、及び流体に接触する壁を表した第２の粒子毎に、速度、位置、圧力に関するデータを記録しておき、接触角を加味して流体の挙動を解析する解析装置が知られている（例えば、特許文献１参照）。【０００６】
ところで、粒子法の計算では、連続体を粒子分布に分けて、連続体の応力による合力を粒子間の相互作用の力の和として表現する。そして、粒子同士が互いに近づこうとする場合には、圧力が大きくなって反発力が発生し、粒子が余りに密集したり疎になったりしない仕組みになっている。【０００７】
相互作用は、粒子の分布から連続場を構成するのに用いられるカーネル関数による重ね合わせにより表現される。例えばＳＰＨ法の場合、流体の運動を解く問題は、例えば、次式（１）〜（４）で表される常微分方程式系に帰着される。式（１）、（２）は、流体の運動量保存則を離散化したものであり、式（３）は、質量保存則を離散化したものである。また、式（４）は、流体の方程式であり、ここでは密度と圧力が線形の関係となる単純な例を挙げているが、式（４）とは異なる状態方程式を採用することも出来る。式中、添え字ｉは、粒子の番号を表している。ｉ番目の粒子の位置ベクトル、速度ベクトル、密度、圧力はそれぞれ#ｒｉ、#ｖｉ、ρｉ、ｐｉである。ここで、「#」は、後に続くアルファベットがベクトルであることを示すものとする。また、ｍｊはｊ番目の粒子の質量を示している。【０００８】
【数１】
2013042234.tif 000003
【０００９】
上式（３）におけるＷには、カーネル関数として３次のスプライン関数などが用いられる。当該関数が次式（５）で表される場合、ｈはカーネル関数が値をもつ範囲を表す。【００１０】
【数２】
2013042234.tif 000004
【００１１】
ｈは粒子間の影響半径であり、初期状態の平均粒子間隔の２倍から３倍程度が良く用いられる。βはカーネル関数の全空間積分量が１になるように調整された値であり、二次元の場合は０．７πｈ２、三次元の場合はπｈ３と決められる。【００１２】
ＳＰＨ法では、式（３）より、粒子同士が近づこうとすると密度が大きくなろうとし、式（４）より、密度が大きくなると圧力が大きくなり、式（２）より、圧力が大きい方から小さい方へと力が働くため、粒子が反発力を得るという仕組みになっている。図１は、接近する粒子ｉ、ｊが互いに反発力ｆを得る様子を示す図である。【００１３】
このように、従来の粒子法による計算手法では、境界に固定した粒子を配置することで、境界に近づこうとする粒子に反発力を発生させ、不透性条件（物体が他の物体の境界を突き抜けないこと）を表現している。流体-構造連成の場合には、構造物の境界粒子が受ける圧力に基づいて、流体と構造物の相互作用を表現している。図２は、流体を構成する粒子Ｄ１が、壁の粒子Ｅ１、Ｅ２、…から反発力ｆを受ける様子を示す図である。
【先行技術文献】
【特許文献】
【００１４】
【特許文献１】特開２００８−１１１６７５号公報
【非特許文献】
【００１５】
【非特許文献１】S.Koshizuka,  A.Nobe  and  Y.Oka,  “Numerical  Analysis  of  Breaking  Waves  Using  The  Moving  Particle  Semi-Implicit  Method”,  International  Journal  for  Numerical  Method  in  Fluids,  26:751-769（1998）
【非特許文献２】J.J.Monaghan,  “Smoothed Particle Hydrodynamics”,  Annu.  Rev.  Astron.  Astrophys.  30:  543-74（1992）
【発明の概要】
【発明が解決しようとする課題】
【００１６】
しかしながら、粒子法を用いた従来の手法では、流体の粒子が構造物の境界を越えない、という境界条件を保障することができない。図３は、粒子法を用いた従来の手法において、境界条件を満たさなくなる場面を示す図である。従来の手法では、壁側の物体が伸展、変形等して粒子間隔が大きくなり、接近する粒子Ｄ１の影響半径ｈの２倍を超えた場合、壁側の物体を構成する粒子から十分な反発力を受けず、境界を突き抜けてしまうという非物理的な計算結果が出力されてしまう。【００１７】
１つの側面では、解析結果が非物理的な挙動を示してしまうことを抑制することが可能な物体運動解析装置等を提供することを目的とする。
【課題を解決するための手段】
【００１８】
上記目的を達成するための本発明の一態様は、第１の物体を表現した複数の粒子を含む第１の粒子群の物理量に関するデータに基づき、前記第１の粒子群に含まれる前記複数の粒子の間を接続する接続要素を生成する生成部と、第２の物体に含まれる第２の粒子の物理量に関するデータを参照し、前記第２の粒子と前記生成された接続要素との最短距離を算出して、前記第１の物体と前記第２の粒子との間に作用する反発力を算出する算出部と、前記算出された反発力を用いて、前記第１の物体及び前記第２の物体の運動を解析する解析部と、を有することを特徴とする物体運動解析装置である。
【発明の効果】
【００１９】
１つの側面では、解析結果が非物理的な挙動を示してしまうことを抑制することが可能な、物体運動解析装置、物体運動解析方法、及び物体運動解析プログラムを提供することができる。
【図面の簡単な説明】
【００２０】
【図１】接近する粒子ｉ、ｊが互いに反発力ｆを得る様子を示す図である。
【図２】流体を構成する粒子Ｄ１が、壁の粒子Ｅ１、Ｅ２、…から反発力ｆを受ける様子を示す図である。
【図３】粒子法を用いた従来の手法において、境界条件を満たさなくなる場面を示す図である。
【図４】本発明の第１実施例に係る物体運動解析装置１のハードウエア構成例である。
【図５】本発明の第１実施例に係る物体運動解析装置１の機能構成例である。
【図６】物体Ｅの外縁に位置する粒子からパネル要素ｉを設定する様子を示す図である。
【図７】Ｄ#ｐ，ｉが０から１の間にある場合の位置関係及び最短距離を示す図である。
【図８】Ｄ#ｐ，ｉが０から１の間に無い場合の位置関係及び最短距離を示す図である。
【図９】ポテンシャル関数ψの性質を示す図である。
【図１０】粒子Ｄ１と物体Ｅとの間に作用する反発力ｆＤ１Ｅを模式的に示す図である。
【図１１】本発明の第２実施例に係る物体運動解析装置２の機能構成例である。
【図１２】物体Ｅの外縁に位置する粒子からパネル要素ｉ（ｊ１２３、ｊ２３４、…）を設定する様子を示す図である。
【図１３】点#ｐから三角形パネル要素ｉを含む平面への射影点を求める様子を示す図である。
【図１４】点#ｐから三角形パネル要素ｉを含む平面への射影点が存在する領域を場合分けした図である。
【図１５】第２実施例の反発力算出部４４により実行される処理の流れを示すフローチャートである。
【図１６】本発明の第３実施例に係る物体運動解析装置３の機能構成例である。
【図１７】物体運動解析装置３が扱う二次元の流体−構造連成問題を模式的に示す図である。
【図１８】ゴム膜を分割して分点（粒子）を定義し、これに基づき分点（粒子）を接続するパネル要素を生成する様子を示す図である。
【図１９】本実施例の物体運動解析装置３により実行される処理の流れを示すフローチャートである。
【図２０】本適用例において説明した手法を用いて実際にシミュレーションを行った際の各過程の状態を示す図である。
【発明を実施するための形態】
【００２１】
以下、本発明を実施するための形態について、添付図面を参照しながら実施例を挙げて説明する。第１実施例では二次元平面上での運動解析について説明し、第２実施例では三次元空間上での運動解析について説明し、第３実施例では、二次元の流体−構造連成問題における運動解析について説明する。
【００２２】
なお、図面及び数式内ではベクトル表記を採用するが、文章中では、「#」が、後に続くアルファベットがベクトルであることを示すものとする。
【実施例】
【００２３】
＜第１実施例＞
以下、本発明の第１実施例に係る物体運動解析装置、物体運動解析方法及び物体運動解析プログラムについて説明する。【００２４】
［ハードウエア構成］
図４は、本発明の第１実施例に係る物体運動解析装置１のハードウエア構成例である。物体運動解析装置１は、例えば、ＣＰＵ（Central  Processing  Unit）１０と、ドライブ装置１２と、補助記憶装置１６と、メモリ装置１８と、インターフェース装置２０と、入力装置２２と、ディスプレイ装置２４と、を備える情報処理装置である。これらの構成要素は、バスやシリアル回線等を介して接続されている。【００２５】
ＣＰＵ１０は、例えば、プログラムカウンタや命令デコーダ、各種演算器、ＬＳＵ（Load  Store  Unit）、汎用レジスタ等を有するプロセッサである。【００２６】
ドライブ装置１２は、記憶媒体１４からプログラムやデータを読み込み可能な装置である。プログラムを記録した記憶媒体１４がドライブ装置１２に装着されると、プログラムが記憶媒体１４からドライブ装置１２を介して補助記憶装置１６にインストールされる。記憶媒体１４は、例えば、ＣＤ−ＲＯＭ（Compact  Disk  Read  Only  Memory）、ＤＶＤ（Digital  Versatile  Disk）、ＵＳＢ（Universal  Serial  Bus）メモリ等の可搬型の記憶媒体である。また、補助記憶装置１６は、例えば、ＨＤＤ（Hard  Disk  Drive）やフラッシュメモリである。【００２７】
プログラムのインストールは、上記のように記憶媒体１４を用いる他、インターフェース装置２０がネットワークを介して他のコンピュータよりダウンロードし、補助記憶装置１６にインストールすることによって行うこともできる。ネットワークは、インターネット、ＬＡＮ（Local  Area  Network）、無線ネットワーク等である。また、プログラムは、情報処理装置の出荷時に、予め補助記憶装置１６やＲＯＭ（Read  Only  Memory）等に格納されていてもよい。【００２８】
このようにしてインストール又は予め格納されたプログラムをＣＰＵ１０が実行することにより、図４に示す態様の情報処理装置が、本実施例の物体運動解析装置１として機能することができる。【００２９】
メモリ装置１８は、例えば、ＲＡＭ（Random  Access  Memory）やＥＥＰＲＯＭ（Electrically  Erasable  and  Programmable  Read  Only  Memory）である。インターフェース装置２０は、上記ネットワークとの接続等を制御する。【００３０】
入力装置２２は、例えば、キーボードやマウス、タッチパッド、タッチパネル、マイク等である。また、ディスプレイ装置２４は、例えば、ＬＣＤ（Liquid  Crystal  Display）やＣＲＴ（Cathode  Ray  Tube）等の表示装置である。物体運動解析装置１は、ディスプレイ装置２４の他、プリンタ、スピーカ等の他の種類の出力装置を備えてもよい。【００３１】
［機能構成］
図５は、本発明の第１実施例に係る物体運動解析装置１の機能構成例である。物体運動解析装置１は、パラメータ入力部３０と、パネル要素生成部３２と、反発力算出部３４と、運動解析部３６と、を備える。これらの機能ブロックは、補助記憶装置１６等に格納されたプログラム・ソフトウエアをＣＰＵ１０が実行することにより機能する。なお、これらの機能ブロックが明確に分離したプログラムによって実現される必要はなく、サブルーチンや関数として他のプログラムによって呼び出されるものであってもよい。また、機能ブロックの一部が、ＩＣ（Integrated  Circuit）やＦＰＧＡ（Field  Programmable  Gate  Array）等のハードウエア手段であっても構わない。【００３２】
パラメータ入力部３０は、解析対象の物体を構成する粒子の座標、速度、結合関係、質量等に関する物理量のデータの入力を受け付け、メモリ装置１８に格納する。これらのデータは、予め補助記憶装置１６に格納されていたものであってもよいし、ユーザが入力装置２２を用いて入力したものであってもよい。また、これらのデータは、記憶媒体１４から補助記憶装置１６にインストールされたものであってもよいし、インターフェース装置２０によりネットワークから取得されたものであってもよい。【００３３】
入力される解析対象の物体の態様は如何なる態様であってもよい。本実施例では、特に、ある平面上に存在する物体Ｄに含まれる粒子Ｄ１と、同一平面上に存在する他の物体Ｅに含まれる複数の粒子（Ｅ１、Ｅ２、Ｅ３、…）と、の間に作用する反発力を算出する手法について説明する。係る反発力以外の物理量（重力や速度、回転速度等）については、他の一般的な物理的手法を用いて算出すればよい。【００３４】
なお、以下では物体Ｄが単一粒子の物体であるかのように説明するが、例えば第３実施例で説明するように、物体Ｄに含まれる複数の粒子が、それぞれ反発力を物体Ｅから受けるものとして運動を解析することも可能である。【００３５】
パネル要素生成部３２は、メモリ装置１８に格納されたデータを読み込み、物体Ｅの外縁（境界線）に位置する粒子（Ｅ１、Ｅ２、Ｅ３、…）を抽出してメモリ装置１８に格納する。そして、パネル要素生成部３２は、抽出した粒子を直線で連結し、粒子間を接続する各線分を接続要素である各パネル要素ｉと定義してメモリ装置１８に格納する。図６は、物体Ｅの外縁に位置する粒子からパネル要素ｉ（ｉ１２、ｉ２３、ｉ３４、…）を設定する様子を示す図である。【００３６】
以下、パネル要素ｉは、始点の座標#ｒｓ，ｉと、終点の座標#ｒｅ，ｉで特定されるものとする。【００３７】
反発力算出部３４は、まず、粒子Ｄ１と、上記設定した各パネル要素との最短距離を算出する。ここでは、粒子Ｄ１の位置を点#ｐとする。反発力算出部３４は、例えば次式（６）〜（８）で表される量を算出し、粒子Ｄ１とパネル要素ｉとの位置関係を定義する。【００３８】
【数３】
2013042234.tif 000005
【００３９】
式（８）におけるＤ#ｐ，ｉが０＜Ｄ#ｐ，ｉ＜１の範囲にある場合は、点#ｐからパネル要素ｉへの垂線が、パネル要素ｉと交わるような位置関係にある。この場合、反発力算出部３４は、パネル要素ｉに落とした垂線の長さを、点#ｐとパネル要素ｉとの最短距離として採用する。一方、反発力算出部３４は、Ｄ#ｐ，ｉがＤ#ｐ，ｉ≦０又は１≦Ｄ#ｐ，ｉの範囲にある場合は、パネル要素ｉの両端である点#ｒｓ，ｉと点#ｒｅ，ｉのうち点#ｐから近い方の点と点#ｐとの距離を、点#ｐとパネル要素ｉとの最短距離として採用する。【００４０】
図７は、Ｄ#ｐ，ｉが０＜Ｄ#ｐ，ｉ＜１の範囲にある場合の位置関係及び最短距離を示す図であり、図８は、Ｄ#ｐ，ｉがＤ#ｐ，ｉ≦０又は１≦Ｄ#ｐ，ｉの範囲にある場合の位置関係及び最短距離を示す図である。【００４１】
点#ｐとパネル要素ｉとの最短距離ｄ#ｐ，２（#ｒｓ，ｉ，#ｒｅ，ｉ）は、次式（９）、（１０）によって求められる。【００４２】
【数４】
2013042234.tif 000006
【００４３】
反発力算出部３４は、各パネル要素（ｉ＝１〜ｎ）について、点#ｐとの最短距離ｄ#ｐ，２（#ｒｓ，ｉ，#ｒｅ，ｉ）を計算する。そして、反発力算出部３４は、次式（１１）により、各パネル要素（ｉ＝１〜ｎ）についての点＃ｐとの最短距離ｄ#ｐ，２（#ｒｓ，ｉ，#ｒｅ，ｉ）のうち最も値の小さい距離ｄ#ｐ，ｍｉｎ，２を抽出し、距離ｄ#ｐ，ｍｉｎ，２に基づき次式（１２）によりポテンシャル関数ψを定義する。図９は、ポテンシャル関数ψの性質を示す図である。図示するようにポテンシャル関数ψは、距離ｄ#ｐ，ｍｉｎ，２がゼロに近づくと値及び傾斜が大きくなり、距離ｄ#ｐ，ｍｉｎ，２がｒ０に近づくとゼロに漸近する。更に、反発力算出部３４は、ポテンシャル関数ψを微分することによって、粒子Ｄ１と物体Ｅとの間に作用する反発力ｆＤ１Ｅを算出する。図９に示すポテンシャル関数ψの性質から、反発力ｆＤ１Ｅは、距離ｄ#ｐ，ｍｉｎ，２が減少してゼロに近づくと値が大きくなり、距離ｄ#ｐ，ｍｉｎ，２が増加してｒ０に近づくとゼロに近づく傾向を有している。図１０は、粒子Ｄ１と物体Ｅとの間に作用する反発力ｆＤ１Ｅを模式的に示す図である。図中、ｉ４５が、粒子Ｄ１との最短距離が最も小さい（すなわち粒子Ｄ１に最も近い）パネル要素である。【００４４】
【数５】
2013042234.tif 000007
【００４５】
運動解析部３６は、反発力算出部３４によって算出された反発力ｆＤ１Ｅ、及びその他の物理量（重力や速度、回転速度等）を用いて、物体Ｄ及びＥの運動を解析する。物体に作用する力に基づく運動解析については、既に種々の手法が公知となっているため、説明を省略する。【００４６】
ここで、従来の手法との比較について説明する。前述したように、粒子法を用いた従来の手法では、壁側の物体が伸展、変形等して粒子間隔が大きくなると、接近する粒子が壁側の物体を構成する粒子から十分な反発力を受けず、境界を突き抜けてしまうという非物理的な計算結果が出力されてしまう場合がある。また、これを避けるために、構造物の境界に多くの粒子を配置すると、計算量が大幅に増加してしまう。【００４７】
この点、本実施例の物体運動解析装置１では、物体Ｄから、物体Ｅの粒子間を接続するパネル要素までの最短距離に基づいて反発力を計算するため、仮に物体Ｅが大きく変形した場合でも、物体が境界を突き抜けるような計算結果は出力されない。また、物体Ｅの外縁（境界線）に多くの粒子を設定しなくとも物体Ｅの外縁をパネル要素で隙間無く覆うことができるため、計算コストの増加を抑制することができる。従って、非物理的な挙動を防止しつつ迅速に結果を出力することができる。【００４８】
以上説明した本実施例の物体運動解析装置、物体運動解析方法及び物体運動解析プログラムによれば、非物理的な挙動を防止しつつ迅速に結果を出力することができる。【００４９】
＜第２実施例＞
以下、本発明の第２実施例に係る物体運動解析装置、物体運動解析方法及び物体運動解析プログラムについて説明する。【００５０】
［ハードウエア構成］については、第１実施例と同様であるため、説明を省略する。【００５１】
［機能構成］
図１１は、本発明の第２実施例に係る物体運動解析装置２の機能構成例である。物体運動解析装置２は、パラメータ入力部４０と、パネル要素生成部４２と、反発力算出部４４と、運動解析部４６と、を備える。これらの機能ブロックは、補助記憶装置１６等に格納されたプログラム・ソフトウエアをＣＰＵ１０が実行することにより機能する。なお、これらの機能ブロックが明確に分離したプログラムによって実現される必要はなく、サブルーチンや関数として他のプログラムによって呼び出されるものであってもよい。また、機能ブロックの一部が、ＩＣ（Integrated  Circuit）やＦＰＧＡ（Field  Programmable  Gate  Array）等のハードウエア手段であっても構わない。【００５２】
パラメータ入力部４０は、物体を構成する粒子の座標、速度、結合関係、質量等に関するデータの入力を受け付け、メモリ装置１８に格納する。本実施例のパラメータ入力部４０が受け付ける物理量のデータは、予め補助記憶装置１６に格納されていたものであってもよいし、ユーザが入力装置２２を用いて入力したものであってもよい。また、これらのデータは、記憶媒体１４から補助記憶装置１６にインストールされたものであってもよい。入力される物体の態様は如何なる態様をもとり得るが、本実施例では、特に、ある空間に存在する物体Ｄに含まれる粒子Ｄ１と、同一空間に存在する他の物体Ｅに含まれる複数の粒子（Ｅ１、Ｅ２、Ｅ３、…）と、の間に作用する反発力を算出する手法について説明する。係る反発力以外の物理量（重力や速度、回転速度等）については、他の一般的な物理的手法を用いて算出すればよい。【００５３】
パネル要素生成部４２は、メモリ装置１８に格納されたデータを読み込み、物体Ｅの外縁（境界面）に位置する粒子（Ｅ１、Ｅ２、Ｅ３、…）を抽出してメモリ装置１８に格納する。そして、パネル要素生成部３２は、抽出した粒子をメッシュ状に連結し、三粒子間を接続する線分によって区画される三角形の面を、接続要素であるパネル要素ｉと定義してメモリ装置１８に格納する。図１２は、物体Ｅの外縁に位置する粒子からパネル要素ｉ（ｉ１２３、ｉ２３４、…）を設定する様子を示す図である。なお、三次元で解析を行う場合のパネル要素ｉは、三角形に限らず、四角形その他の形状を有してもよい。【００５４】
以下、パネル要素ｉは、三点の座標#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉで特定されるものとする。【００５５】
反発力算出部４４は、まず、粒子Ｄ１と、上記設定された各パネル要素との最短距離を算出する。以下、第１実施例と同様、粒子Ｄ１の位置を点#ｐとする。反発力算出部４４は、以下の手順で、点#ｐと、三点で構成される閉領域である三角形パネル要素ｉとの最短距離を算出する。【００５６】
まず、反発力算出部４４は、点#ｐから三角形パネル要素ｉを含む平面への射影点#ｈｉを求める。図１３は、点#ｐから三角形パネル要素ｉを含む平面への射影点#ｈｉを求める様子を示す図である。【００５７】
反発力算出部４４は、射影点#ｈｉが三角形パネル要素ｉに含まれる場合は、点#ｐから射影点#ｈｉまでの長さを最短距離とする。【００５８】
一方、反発力算出部４４は、射影点#ｈｉが三角形パネル要素ｉに含まれない場合は、点#ｐから、三角形パネル要素ｉの外縁線において射影点#ｈｉに最も近い点までの距離を最短距離とする。【００５９】
以下、上記の手順を数式により表す。まず、点#ｐから三角形パネル要素ｉを含む平面への射影点#ｈｉは、次式（１３）により表される。【００６０】
【数６】
2013042234.tif 000008
【００６１】
ここで、ｄ#ｒ２，ｉ＝ｒ２，ｉ−ｒ１，ｉ、ｄ#ｒ３，ｉ＝ｒ３，ｉ−ｒ１，ｉとすると、#ｐ−＃ｈｉは、ｄ#ｒ２，ｉとｄ#ｒ３，ｉに対して垂直であるから、次式（１４）、（１５）が成立する。この両式を満たすようにパラメータｓ２、ｓ３を決定することで、射影点#ｈｉを求めることができる。【００６２】
【数７】
2013042234.tif 000009
【００６３】
反発力算出部４４は、射影点#ｈｉを求めると、以下のように場合分けして、各点#ｐと三角形パネル要素ｉとの最短距離を算出する。【００６４】
（１）ｓ２≧０，ｓ３≧０，且つｓ２＋ｓ３≧１の場合は、射影点#ｈｉが三角形パネル要素ｉに含まれるため、点#ｐと各三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）は、次式（１６）で表される。【００６５】
【数８】
2013042234.tif 000010
【００６６】
（２）ｓ２＜０，ｓ３＜０，且つｓ２＋ｓ３≧１の場合は、射影点#ｈｉが図１４における領域Ｑ１に存在することになる。図１４は、粒子Ｅ１の位置から各三角形パネル要素ｉを含む平面への射影点が存在する領域を場合分けした図である。この場合、点#ｐと各三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）を求める問題は、点#ｐと線分の最短距離を求める問題となる。従って、反発力算出部４４は、第１実施例の手法を援用し、#ｒｓ，ｉ＝#ｒ２，ｉ、#ｒｅ，ｉ＝#ｒ３，ｉと代入して求めた次式（１７）の値を、点#ｐと接続要素である三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）とする。【００６７】
【数９】
2013042234.tif 000011
【００６８】
（３）ｓ２≧０，ｓ３＜０，且つｓ２＋ｓ３≧１の場合は、射影点#ｈｉが図１４における領域Ｑ２に存在することになる。この場合も、点#ｐと各三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）を求める問題は、点#ｐと線分の最短距離を求める問題となる。従って、反発力算出部４４は、第１実施例の手法を援用し、#ｒｓ，ｉ＝#ｒ１，ｉ、#ｒｅ，ｉ＝#ｒ２，ｉと代入して求めた次式（１８）の値を、点#ｐと接続要素である三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）とする。【００６９】
【数１０】
2013042234.tif 000012
【００７０】
（４）ｓ２＜０，ｓ３≧０，且つｓ２＋ｓ３≧１の場合は、射影点#ｈｉが図１４における領域Ｑ３に存在することになる。この場合も、点#ｐと各三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）を求める問題は、点#ｐと線分の最短距離を求める問題となる。従って、反発力算出部４４は、第１実施例の手法を援用し、#ｒｓ，ｉ＝#ｒ３，ｉ、#ｒｅ，ｉ＝#ｒ１，ｉと代入して求めた次式（１９）の値を、点#ｐと接続要素である三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）とする。【００７１】
【数１１】
2013042234.tif 000013
【００７２】
（５）上記（１）〜（４）以外の場合は、射影点#ｈｉと各三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）は、射影点#ｈｉと各三角形パネル要素ｉのいずれかの頂点との距離となる。従って、反発力算出部４４は、次式（２０）に基づき点#ｐと各三角形パネル要素ｉとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）を算出する。【００７３】
【数１２】
2013042234.tif 000014
【００７４】
反発力算出部４４は、以上の手順によって、各三角形パネル要素ｉについて、点#ｐとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）を算出する。そして、次式（２１）により、各三角形パネル要素ｉについての点＃ｐとの最短距離ｄ#ｐ，３（#ｒ１，ｉ，#ｒ２，ｉ，#ｒ３，ｉ）のうち、最も値の小さい距離ｄ#ｐ，ｍｉｎ，３を算出する。更に、反発力算出部４４は、距離ｄ#ｐ，ｍｉｎ，３に基づき、次式（２２）によりポテンシャル関数ψを定義し、ポテンシャル関数ψを微分することによって、粒子Ｄ１と物体Ｅとの間に作用する反発力ｆＤ１Ｅを算出する。【００７５】
【数１３】
2013042234.tif 000015
【００７６】
運動解析部４６は、反発力算出部４４によって算出された反発力ｆＤ１Ｅ、及びその他の物理量（重力や速度、回転速度等）を考慮し、物体Ｄ及びＥの運動を解析する。物体に作用する力に基づく運動解析については、既に種々の手法が公知となっているため、説明を省略する。【００７７】
以下、上記説明した反発力算出部４４の処理を、フローチャートの形式で説明する。図１５は、第２実施例の反発力算出部４４により実行される処理の流れを示すフローチャートである。本フローでは、物体Ｄが複数の粒子（Ｄ１、Ｄ２、Ｄ３、…）を含むものとする。【００７８】
反発力算出部４４は、Ｓ１００〜Ｓ１１２の処理を、物体Ｄを構成する各粒子Ｄｊ（ｊ＝１〜ｋ）について実行する。また、反発力算出部４４は、Ｓ１００〜Ｓ１０８の処理を、各三角形パネル要素ｉ（ｉ＝１〜ｎ）について実行する。【００７９】
反発力算出部４４は、粒子Ｄｊから三角形パネル要素ｉを含む平面への射影点＃ｈｉを求める（Ｓ１００）。【００８０】
次に、反発力算出部４４は、射影点が三角形パネル要素ｉに含まれるか否かを判定する（Ｓ１０２）。【００８１】
反発力算出部４４は、射影点が三角形パネル要素ｉに含まれる場合には、粒子Ｄｊから射影点までの距離を、粒子Ｄｊから三角形パネル要素ｉまでの最短距離とする（Ｓ１０４）。【００８２】
反発力算出部４４は、射影点が三角形パネル要素ｉに含まれない場合には、三角形パネル要素ｉの外縁線（頂点を含む）上で、粒子Ｄｊから最も近い点を抽出する（Ｓ１０６）。次に、反発力算出部４４は、粒子Ｄｊから抽出した点までの距離を、粒子Ｄｊから三角形パネル要素ｉまでの最短距離とする（Ｓ１０８）。【００８３】
反発力算出部４４は、Ｓ１００〜Ｓ１０８の処理を全ての三角形パネル要素ｉについて実行すると、各三角形パネル要素ｉについて算出した最短距離の中で、最も値が小さいものを抽出する（Ｓ１１０）。【００８４】
そして、反発力算出部４４は、抽出した値を用いて、粒子Ｄｊと物体Ｅと間の反発力を表すポテンシャル関数を算出する（Ｓ１１２）。【００８５】
このような処理によって奏される効果は、第１実施例と同様である。従って、従来技術との比較等については、第１実施例の説明を三次元に拡張すればよいため、説明を省略する。【００８６】
以上説明した本実施例の物体運動解析装置、物体運動解析方法、及び物体運動解析プログラムによれば、非物理的な挙動を防止しつつ迅速に結果を出力することができる。【００８７】
＜第３実施例＞
以下、本発明の第３実施例に係る物体運動解析装置、物体運動解析方法、及び物体運動解析プログラムについて説明する。第３実施例の物体運動解析装置３は、第１実施例の物体運動解析装置１が有する機能を基本とし、二次元の流体−構造連成問題に適用可能となるように機能を拡張したものである。【００８８】
ハードウエア構成については、第１実施例や第２実施例と同様であるため、説明を省略する。【００８９】
［機能構成］
図１６は、本発明の第３実施例に係る物体運動解析装置３の機能構成例である。物体運動解析装置３は、パラメータ入力部５０と、粒子点設定部５１と、パネル要素生成部５２と、反発力算出部５４と、運動解析部５６と、を備える。これらの機能ブロックは、補助記憶装置１６等に格納されたプログラム・ソフトウエアをＣＰＵ１０が実行することにより機能する。なお、これらの機能ブロックが明確に分離したプログラムによって実現される必要はなく、サブルーチンや関数として他のプログラムによって呼び出されるものであってもよい。また、機能ブロックの一部が、ＩＣ（Integrated  Circuit）やＦＰＧＡ（Field  Programmable  Gate  Array）等のハードウエア手段であっても構わない。なお、各機能ブロックの機能については、後述する図１９のフローに則して説明する。【００９０】
物体運動解析装置３は、二次元のゴム膜（第１実施例における物体Ｅの一例）と流体（第１実施例における物体Ｄの一例）が相互作用する支配方程式系の問題を取り扱うことができる。図１７は、物体運動解析装置３が扱う二次元の流体−構造連成問題を模式的に示す図である。図示するように、物体運動解析装置３は、ゴム膜Ｅの上に純水等の流体Ｄが置かれた系における、ゴム膜Ｅや流体Ｄの運動を解析する。以下、この系における高さをｙ座標で、横方向の位置をｘ座標で表すものとする。【００９１】
まず、物体運動解析装置３が扱う支配方程式を、次式（２３）〜（２６）で定義する。式中、ｕ（ｘ，ｔ）はゴム膜を構成する各粒子の時間ｔにおける位置ｘに対する高さを表す関数であり、ρ（ｘ，ｙ，ｔ）、ｐ（ｘ，ｙ，ｔ）、#ｖ（ｘ，ｙ，ｔ）はそれぞれ、流体の密度場、圧力場、速度場である。#Ｆは外力であり、ここでは重力を表すものとする。式（２３）はゴム膜の運動量保存則であり、式（２４）は流体の質量保存則であり、式（２５）は流体の運動量保存則であり、式（２６）は流体の状態方程式である。【００９２】
【数１４】
2013042234.tif 000016
【００９３】
また、Ｐｕは、流体がゴム膜を押す力を表す。この項は以下のように導出される。（ｘ−ｄｘ／２）と（ｘ＋ｄｘ／２）の間にある膜の微小要素に作用する力は、−ｐ（ｘ，ｕ（ｘ））ｄ#Ｓ・#ｙ＾で表される。ここで、ｄ#Ｓは上記微小要素の面積を有し、向きが外向き法線方向の面積ベクトルであり、#ｙ＾はｙ方向の単位ベクトルである。また、ｄ#Ｓ＝（−ｕｘ，１）ｄｘであるため、Ｐｕは、単位長さあたりの力として表現すると、次式（２７）で表される。【００９４】
Ｐｕ＝−ｐ（ｘ，ｕ（ｘ））  …（２７）。【００９５】
ここで、∂ｕ／∂ｔ＝ηとし、ｕを、空間についてはｄｘでｍ等分すると共に、時間についてはｄｔで一様に分割し、ｕ（ｉ×ｄｘ，ｎ×ｄｔ）＝ｕｉｎという表記を用いると、上式（２３）を、次式（２８）〜（３０）のように離散化することができる。式中、#ｕ*＝（ｕ０*，ｕ１*，…，ｕｉ*，…，ｕｍ*）は、各分点における値ｕｉ*を並べたベクトルである。また、ψｂは、流体を構成する任意の粒子ｂと、パネル要素との距離から得られるポテンシャル関数である。また、上付添え字「*」は、時間に関する中間値を表している。なお、式（２９）の第２項は、張力係数γが乗算されたゴム膜の内部力を示している。【００９６】
【数１５】
2013042234.tif 000017
【００９７】
図１８は、ゴム膜を分割して分点（粒子）を定義し、これに基づき分点（粒子）を接続するパネル要素を生成する様子を示す図である。図１８に示すように、高さがｕ（ｘ，ｔ）で表されるゴム膜は、例えばｄｘ幅を持つようにｍ等分され、分割された切片の中心座標が仮想的な粒子の座標（ｕ１、ｕ２、…ｕｉ、ｕｉ＋１、…）と認識される。そして、分割認識された仮想的な粒子を連結する各線分が、接続要素である各パネル要素ｉとして認識及び処理される。
点（ｉ×ｄｘ，ｕｉ*）と点（（ｉ＋１）×ｄｘ，ｕｉ＋１*）で作られる各線分要素をｉ番目のパネル要素とすると、ゴム膜は、ｍ−１個のパネル要素で構成される。
【数１６】
2013042234.tif 000018
として実施例１の手順でポテンシャル関数ψｂを定義すると、ψｂは#ｒｂ*と#ｕ*のみの関数となるため、ポテンシャル関数ψｂをψｂ（#ｒｂ*，#ｕ*）と表記した。【００９８】
ここで、ポテンシャル関数ψｂの微分∂ψｂ（#ｒｂ*，#ｕ*）／∂ｕｉ*を解析的に求めるのは困難なため、本実施例では、次式（３１）に示す差分近似により求める。式中、ｕｉ，＋*＝（ｕ０*，ｕ１*，…，ｕｉ*＋δ，…，ｕｍ*）、ｕｉ，−*＝（ｕ０*，ｕ１*，…，ｕｉ*−δ，…，ｕｍ*）である。また、δは、例えば１０−８等の定数である。次式（３１）は、パネル要素が流体から受ける反発力を示している。【００９９】
【数１７】
2013042234.tif 000019
【０１００】
一方、流体の方程式（２４）〜（２６）は、ＳＰＨ法に基づき、次式（３２）〜（３５）により離散化される。式中、βは定数である。また、Ｌ２（#ｒａ*）は二次元の再標準化行列であり、次式（３６）で表される。【０１０１】
【数１８】
2013042234.tif 000020
【０１０２】
ここでは、
【数１９】
2013042234.tif 000021
が成立する。ｐａｂｎ＋１／２、ｖａｂｎ＋１／２は、流体を構成する任意の粒子ａｂ間での一次元リーマン問題を解いて求めた時空間の中間値である。ｐａｂｎ＋１／２、ｖａｂｎ＋１／２は、具体的には、以下のように決定される。【０１０３】
まず流体に含まれる粒子ａ、ｂに関して、次式（３７）〜（４０）で表される特性関数を定める。【０１０４】
【数２０】
2013042234.tif 000022
【０１０５】
更に、勾配を次式（４１）〜（４６）により計算する。式中、速度の上付き添え字は、各座標の成分を表している。【０１０６】
【数２１】
2013042234.tif 000023
【０１０７】
これらを用いて、ｐａｂｎ＋１／２、ｖａｂｎ＋１／２を次式（４７）〜（５１）により決定することができ、上式（３３）に示すポテンシャル関数ψａの微分を、次式（５２）の差分近似によって求めることができる。【０１０８】
【数２２】
2013042234.tif 000024
【０１０９】
ここでは、
【数２３】
2013042234.tif 000025
が成立する。次式（５２）は、流体がパネル要素から受ける反発力を示している。【０１１０】
【数２４】
2013042234.tif 000026
【０１１１】
以上のような計算により、第１実施例で説明した運動解析原理を、二次元の流体−構造連成問題に適用することができる。また、同様に、第２実施例で説明した運動解析原理についても、流体や構造体としてデータが与えられたモデルに対して適用することが可能である。【０１１２】
図１９は、本実施例の物体運動解析装置３により実行される処理の流れを示すフローチャートである。なお、本フローチャートは、構造物と流体を処理対象とするが、構造物同士、流体同士等についても適用可能である。【０１１３】
まず、パラメータ入力部５０が、構造物（例えばゴム膜）及び流体についてのデータを取得する（Ｓ２００）。パラメータ入力部５０により取得されるデータは、構造物の座標、速度、質量、及び結合関係（張力係数等を含む）、並びに流体を構成する粒子の座標、速度、質量等である。【０１１４】
次に、粒子点設定部５１が、上記のように構造物を分割して仮想的な複数の粒子を設定する（Ｓ２０２）。本ステップの処理の概要は、図１８に則して説明した通りである。このように、粒子点設定部５１は、構造物の特に外縁（境界線）を所望の間隔で分割し、仮想的な粒子を設定する。これによって、入力データが粒子として与えられていない場合であっても、パネル要素生成部５２にパネル要素を生成させることができる。この結果、入力データが粒子データとして与えられていない場合でも、第１実施例や第２実施例と同様に、物体が境界を突き抜けるような非物理的な挙動を防止しつつ迅速に結果を出力することができる。【０１１５】
次に、パネル要素生成部５２が、粒子点設定部５１が設定した仮想的な粒子を接続してパネル要素を生成する（Ｓ２０４）。【０１１６】
次に、運動解析部５６が、Ｓ２００において入力された座標及び速度、並びにＳ２０２及びＳ２０４において生成されたパネル要素のデータを用いて、時刻をｄｔ／２進め、構造物及び流体の粒子及びパネル要素等の座標を更新する（Ｓ２０６）。【０１１７】
次に、反発力算出部５４が、Ｓ２０８〜Ｓ２１０の処理を、流体を構成する各粒子（以下、流体粒子と称する）について実行する。【０１１８】
まず、流体粒子を順に一つ取り出し、現時点における各パネル要素との位置関係に基づきポテンシャル関数ψを算出する（Ｓ２０８）。そして、ポテンシャル関数の微分を計算して相互作用項（流体粒子とパネル要素との間に作用する反発力）を算出する（Ｓ２１０）。【０１１９】
運動解析部５６は、各流体粒子とパネル要素との間に作用する反発力がＳ２０８〜Ｓ２１０で算出されると、算出された反発力を構造物の物理モデルに適用して構造物の加速度を算出し（Ｓ２１２）、算出した加速度を用いて構造物の速度を更新する（Ｓ２１４）。【０１２０】
次に、運動解析部５６は、上記算出された反発力を流体の物理モデルに適用して流体粒子の加速度を算出し（Ｓ２１６）、算出した加速度を用いて流体粒子の速度を更新する（Ｓ２１８）。【０１２１】
そして、運動解析部５６は、上記更新した構造物の速度、及び流体粒子の速度を用いて、時刻をｄｔ進め、構造物及び流体粒子の座標を更新する（Ｓ２２０）。【０１２２】
更に、運動解析部５６は、上記更新した構造物及び流体粒子の座標を、ディスプレイ装置２４等を用いて出力し（Ｓ２２２）、Ｓ２０８に戻る。なお、物体運動解析装置３は、Ｓ２０８〜Ｓ２２２のループを所定回数繰り返すと本フローを終了する等と、予め設定されていてもよい。【０１２３】
図２０は、本適用例において説明した手法を用いて実際にシミュレーションを行った際の各過程の状態を示す図である。物体運動解析装置１に与えたパラメータは、Ｌ＝0.292［ｍ］、ｃ＝100［ｍ／ｓ］、Ｍ＝80、ｄｘ＝2Ｌ／Ｍ［ｍ］、ｄｔ＝ｄｘ／（5ｃ）［ｓ］、σ＝10［ｋｇ／ｍ］、α＝3［ｋｇｍ］、γ＝250［Ｎ／ｓ］である。【０１２４】
また、係るシミュレーションでは、流体の初期密度はρａ０＝1000［ｋｇ／ｍ２］で一定にし、流体の粒子の質量は全てρａ０／（ｄｘ）２とし、ｈ＝2ｄｘとし、粒子数を2054個とした。【０１２５】
図示するように、ゴム膜Ｅが大きく伸びた場合でも、流体Ｄがゴム膜Ｅを突き抜けるといった非物理的挙動が起きずに計算を行うことが出来ることが示された。【０１２６】
以上説明した本実施例の物体運動解析装置、物体運動解析方法、及び物体運動解析プログラムによれば、非物理的な挙動を防止しつつ迅速に結果を出力することができる。【０１２７】
また、構造物を分割して仮想的な複数の粒子を設定する粒子点設定部５１を備えることにより、入力データが粒子として与えられていない場合であっても、パネル要素生成部５２にパネル要素を生成させることができる。この結果、入力データが粒子として与えられていない場合でも、第１実施例や第２実施例と同様に、物体が境界を突き抜けるような非物理的な挙動を防止しつつ迅速に結果を出力することができる。【０１２８】
以上、本発明を実施するための最良の形態について実施例を用いて説明したが、本発明はこうした実施例に何等限定されるものではなく、本発明の要旨を逸脱しない範囲内において種々の変形及び置換を加えることができる。【０１２９】
例えば、第１〜第３実施例において、粒子とパネル要素の最短距離のうち最も値の小さい距離を抽出し、これに基づきポテンシャル関数ψを定義するものとした。すなわち、粒子に最も近いパネル要素のみをポテンシャル関数ψの定義に反映させているが、これに限らず、ある粒子と全てのパネル要素についてポテンシャル関数ψを定義して反発力を算出してもよい。
【産業上の利用可能性】
【０１３０】
本発明は、種々の構造物や弾性体、流体の製造産業、製造支援産業、並びにコンピュータソフトウエア産業等に広く利用可能である。
【符号の説明】
【０１３１】
１、２、３        物体運動解析装置
３０、４０、５０  パラメータ入力部
３２、４２、５２  パネル要素生成部
３４、４４、５４  反発力算出部
３６、４６、５６  運動解析部
５１              粒子点設定部
【図１】
2013042234.tif 000027
【図２】
2013042234.tif 000028
【図３】
2013042234.tif 000029
【図４】
2013042234.tif 000030
【図５】
2013042234.tif 000031
【図６】
2013042234.tif 000032
【図７】
2013042234.tif 000033
【図８】
2013042234.tif 000034
【図９】
2013042234.tif 000035
【図１０】
2013042234.tif 000036
【図１１】
2013042234.tif 000037
【図１２】
2013042234.tif 000038
【図１３】
2013042234.tif 000039
【図１４】
2013042234.tif 000040
【図１５】
2013042234.tif 000041
【図１６】
2013042234.tif 000042
【図１７】
2013042234.tif 000043
【図１８】
2013042234.tif 000044
【図１９】
2013042234.tif 000045
【図２０】
2013042234.tif 000046
【国際調査報告】
2013042234.tif 000047
2013042234.tif 000048
2013042234.tif 000049
2013042234.tif 000050