JP2007114866A

JP2007114866A - 接触判定方式

Info

Publication number: JP2007114866A
Application number: JP2005303224A
Authority: JP
Inventors: Toshihiko Sugimoto; 俊彦杉本
Original assignee: Canon Inc
Current assignee: Canon Inc
Priority date: 2005-10-18
Filing date: 2005-10-18
Publication date: 2007-05-10

Abstract

【課題】格子を用いる数値シミュレーションにおいて、精度保証付きの近似距離を用いて複数の物体の接触判定を安価に行う。
【解決手段】二次元または三次元内の格子上に定義された第一の物体、第二の物体の物理量を数値的に模擬する際に必要となる接触判定、及び接触力を算出する方法において、格子間隔でパラメータ化された近似距離を計算する近似距離計算手段を有し、第一の物体、第二の物体夫々の表面上に取った有限個の代表点同士を第一の物体、第二の物体夫々について幾何学要素でつなぎ合わせた近似図形により第一の物体、第二の物体夫々の表面形状を近似し、第二の物体が定義された格子の格子点との近似距離計算手段を用いて計算した符号付き近似距離を格子点において計算することにより所望の領域で符号付き近似距離関数を定義し、近似図形を構成する各代表点における符号付き近似距離関数の値により、接触判定、及び接触力算出を行う。
【選択図】図１

Description

本発明は、接触を伴う構造物の物理量を数値的に模擬及び解析する方法に関する。

構造物の応力、変位などの物理量を数値的に模擬および解析を行う方法は、有限要素法に代表される格子を用いた方法が一般的である。格子は、点、線などの幾何学要素により空間を有限な領域に分割することで、空間の離散化を行うために用いられる。空間分割された格子を構成する該有限領域を以降、胞と呼ぶ。一般的に胞は３次元では多面体、２次元では多角形により構成される場合が多い。この場合、有限要素法においては胞を要素、胞の頂点を節点と呼ぶ。さらに該格子上に取った代表点において、諸々の物理量を離散化して定義する。有限要素法では、節点を代表点として、節点上に諸々の物理量を定義し、要素内の物理量は要素を構成する節点の値から補間を行うことにより算出する。

有限要素法において構造物の解析を行う際、解析対象の物体形状に沿って格子を定義し、物体の移動、変形に合わせて格子も移動、変形させる方法（以降、移動格子法と呼ぶ）が構造物の形状表現精度の面から有利であるなどの理由から広く用いられている。

また、構造物の解析において構造物同士の接触を伴う場合、接触を数値的に扱う方法として、ペナルティ法が広く用いられている。ペナルティ法は接触する物体同士の侵入を許容し、侵入量に応じた力を接触力として物体に課すことで接触を数値的に模擬する方法である。ペナルティ法は接触を外力の境界条件に置き換えて扱うことができるため、その取り扱いが容易である利点がある反面、物体の侵入を許容する手法であるため、接触の表現精度が低下するという欠点もある。ペナルティ法を用いて接触の表現精度を求める場合、小さな侵入量で大きな接触力を与えるようにすることで、該侵入量を小さくすることが可能であるが、この場合、局所的に大きな力のつり合いを計算することになり、過剰に接触の表現精度を求めると、数値計算上の不安定性を招くことになってしまう。そのため、ペナルティ法における接触力の与え方は解析者が求める表現精度と数値計算の安定性を考慮して適切に選択する必要がある。

以上に述べたような理由から、接触を伴う構造物の応力、変位などの物理量を数値的に模擬および解析を行う方法は上記の移動格子法とペナルティ法を組み合わせて用いる方法が広く用いられている。しかしながら、粉体などの多数の物体を該組み合わせ方法により模擬、解析を行う場合、物体同士の接触判定の際、格子同士の幾何学計算により接触判定を行うため、処理が煩雑になり、計算コストも増加するため、このような多数の物体の解析には適した方法ではない。

このような場合、空間に固定された格子を用いた方法（以降、固定格子法と呼ぶ）を用いるのが効果的な方法の一つである。前記の移動格子法では格子形状が近似物体を表現していたが、固定格子法において物体の形状を表現、追跡するために広く用いられているのが符号付き距離関数を用いる方法である。符号付き距離関数は空間上のある点から物体界面までの距離に、その点が物体の内外のどちらの領域に属するかによって正負の符号付けをしたものである。固定格子法ではこの符号付き距離関数を格子上の代表点において定義し、物体の形状を表現する。この方法によると、多数の物体あるいは複雑な形状を持った物体などを扱う場合などでも、一度、全物体に対しての符号付き距離関数を格子の全代表点で求めておけば、符号付き距離関数の符号を調べるだけで容易に接触判定が可能となるなどの利点がある。

ただし、このような距離を計算する手段は、一般に高価であることが知られている。

例えば、２次元平面内に存在するＮ本の線分によって定義された区分的線分曲線Γにおいて、ある点ｐからの区分的線分曲線への距離d（p,Γ）を計算するためには、距離が最小値であることを利用する。つまり、Γを構成する線分(a[i], a[i+1])とｐとの距離を d(p,(a[i], a[i+1]))とすると、
d(p, Γ) = min d (p,(a[i], a[i+1]))
全てのiについて
として、各線分への距離の最小値をもって、d（p,Γ）とする。但しここで、pとa[i],はそれぞれ２次元位置ベクトルでp=(px,py)であり、a[i]＝（a[i]x, a[i]y）としている。

他方、d(p,(a[i], a[i+1]))との距離は、図９のように求まる。線分(a[i], a[i+1])から決まる直線1001への垂線Ｑを求めて、その足Ｑが該線分内であれば、その点まで距離をｐと線分(a[i], a[i+1])とし、線分の外であれば、両端の頂点の近い方の距離を線分への距離とする。垂線の足の決定方法は、(p - λa[i] -(1 -λ)a[i+1], n[i]) =0とする方程式から、λを決定して求まる。

空間が３次元で余次元１の曲面の場合も同様に求まる。その場合も数値計算の場合の多くは物体の界面が区分的線分の代わりに、区分的三角形によるＰＬ（Piecewised Linear）図形によって与えられる。この場合は個々の単体への距離の計算は２次元の場合より、より複雑になるが、計算の手順は上記と同様であり自明である。但し、計算コストは非常に高くなる。

よって符号付き距離関数を利用した接触判定は一般に計算コストが高く、実用に向かない。

従来例としては、例えば特許文献１と特許文献２をあげることが出来る。
特開２００２−５６０３７号公報特開２００２−３５１８５７号公報

接触を伴う多数の物体の模擬、解析方法には固定格子法において符号付き距離関数を用いるのが効果的ではあるが、従来、符号付き距離関数の精度保証なしに、必要以上の過剰精度で計算されているという問題があった。そのため、上記に記したように、計算コストが高いという問題もあった。例えば物質の形状を追跡するような数値計算を、固定格子法を用いて、場の問題として扱う際、該格子のすべての代表点において値が定義されなければならない。このような空間固定の格子を用いて物質の形状を扱う数値計算においては、幾何学量も含めたすべての物理量は格子解像度できまるものとして利用される。つまり、数値計算に利用する立場では、物体の界面からの距離ではなく、格子解像度の数十分の一の解像度による近似距離で十分である。

本発明ではこのような観点から実際の距離ではなく、必要な精度で距離を近似した近似距離関数を使用することにより、所望の精度により、高速に接触を伴う構造物の物理量を数値的に模擬及び解析を行う方法を提供することにある。

本発明においては、二次元または三次元内の格子上に定義された第一の物体、第二の物体の物理量を数値的に模擬する際に必要となる接触判定、及び接触力を算出する方法において、格子間隔でパラメータ化された近似距離を計算する近似距離計算手段を有し、第一の物体、第二の物体それぞれの表面上に取った有限個の代表点同士を第一の物体、第二の物体それぞれについて幾何学要素でつなぎ合わせた近似図形により第一の物体、第二の物体それぞれの表面形状を近似し、該第一の物体、第二の物体と該第一の物体、第二の物体が定義された格子の格子点との該近似距離計算手段を用いて計算した符号付き近似距離を該格子点において計算することにより所望の領域で符号付き近似距離関数を定義し、該近似図形を構成する各代表点における該符号付き近似距離関数の値により、接触判定、及び接触力算出を行う。

本発明によれば、接触を伴う構造体の物理量を数値的に模擬、解析する方法において、格子の格子区間に応じた精度の精度保証付きの高速近似距離計算方法を使用した、物体同士の接触判定方法が提供される。

これにより、多数の物体を模擬、解析対象とする場合などに、特に効率的に十分な精度でかつ高速に接触判定を行うことができる。

（実施例１）
本発明の接触判定方式により接触を伴う場合の構造物の解析を行うプログラムは図１１に示す計算機内のハードディスク（１１７）および、フレキシブルディスク（１１０）に収められている。本実施例では図１に示したフローチャートに従って本発明の接触判定方式により構造物の接触を伴う解析を計算処理装置上で行った例である。

ここで計算処理装置とは図１１に示すようにCPU（１１６）やメモリ（１１８）等からなる演算処理部、ハードディスク（１１７）やフレキシブルディスク１１０からなる記憶装置部、キーボード１１５やマウス１１４等からなる入力部そしてディスプレイ（１１２）等からなる出力表示部などで構成されているシステムである。各演算部等は、本体１１１に収まっている。

計算処理装置はメモリ上に図１に示したフローチャートに沿って設計、コード化されたプログラムをロードし、必要な計算領域を確保し、適当な方法で入力された画像やパラメータに基づき所定の演算処理を行い、その結果得られた領域上での物理量の値をハードディスク等に書き込んで記憶させたり、ディスプレイ上に表示させる。

本発明を実施する一例のアルゴリズムのフローチャートを図１に示す。

なお、本実施例では、図２に示すように空間に構造物Ａ、Bがあり、Ａ、Ｂが接触するときの解析を計算処理装置を用いて行うものとする。

図１のステップＳ１において一連の処理が開始される。ステップＳ２では物体Ａ、Ｂを包含する矩形領域に直交固定格子Ｅを定義する。本実施例において、処理を簡単にするため、均一の直交格子を使用したが、任意の空間固定の格子を使用できる。ステップＳ３ではＡ、Ｂを含むそれぞれの領域に格子Ｌａ，Ｌｂを定義する。Ｌａ，ＬｂはＡ，Ｂの変形、移動を計算するための格子であり、所望の格子を用意すればよい。ステップＳ４ではＡ、Ｂの界面がそれぞれＬａ、Ｌｂと交差する点列Ｔａ、Ｔｂを取得する。この点列Ｔａ，ＴｂによりＡ，Ｂの形状のサンプリングが行われる。ステップＳ５ではＡ，Ｂの接触判定を行うため、Ｔａ、ＴｂをＥ上へ写像する。ステップＳ６ではＴａ、Ｔｂをそれぞれ線分によりつなぎ合わせ、折れ線図形Ｓａ、Ｓｂを定義する。（図１５，図４）ステップＳ７ではＳａ、ＳｂがＥと交差する点列Ｔａ'、Ｔｂ'を取得する。これにより格子サイズに応じた物体形状のサンプリングがなされる。ステップＳ８ではＴａ'、Ｔｂ'それぞれの点列を線分によりつなぎ合わせた折れ線図形Ｓａ'、Ｓｂ'を定義する。（図１６〜１７，図５〜６）ステップＳ９ではＥ上のすべての格子点においてＳａ'、Ｓｂ'の以下で述べる近似距離を利用した符号付き近似距離関数ｆａ、ｆｂを定義する。また同時に格子点を始点、始点における近似距離の定義に採用した始点とペアを成すサンプリング点を終点とするベクトル（以降、近似距離方向ベクトルと呼ぶ。）の定義も行う。（図３，図７）なお、本実施例では領域の内側を負、外側を正として近似距離関数の符号を定めている。

次に、ステップＳ１０ではＳａの各頂点において、ｆｂの値と近似距離方向ベクトルの値、また、Ｓｂの各頂点において、ｆａの値と近似距離方向ベクトルの値を補間する。ステップＳ１１ではＳａ、Ｓｂの各頂点においてのｆｂ、ｆａの値が負であれば接触している頂点（以降、接触頂点と呼ぶ。）と判定される。ここで、接触頂点の近似距離方向ベクトルの大きさを侵入量とする。（図８，図１２）ステップＳ１２では、Ｓａ、Ｓｂの各接触頂点において、近似距離方向ベクトルの向きに侵入量に応じた接触応力をペナルティ法により定義する。ここで、接触応力を定める関数として、σ(x)＝ｋｘ（ｋはペナルテイ係数、ｘは侵入量）なる関数を使用しているが、侵入量に対して単調増加な適当な関数を与えてもよい。また、ペナルティ係数ｋの値として、Ａ，Ｂのヤング率のうち、大きな方の値の１０倍の値を設定している。ステップＳ１３ではＳａ、Ｓｂの各辺の中点において、辺の両端点の値から接触応力を補間し、辺の長さを掛け合わせることで、辺に課される接触力を計算する。ここでは、辺の中点の値を用いて辺が受ける接触力を算出したが、所望の方法により接触応力を辺に沿って線積分を行い、各辺に課される接触力を算出してもよい。ステップＳ１４では、接触力を外力として、所望の方法により、Ａ、Ｂについての力学的な方程式を解くことで、Ａ、Ｂの変形、移動を求める。最後にステップＳ１５において、一連の処理を終了する。

以下において、本実施例のステップＳ９において近似距離及び近似距離方向ベクトルを求める方法を以下で説明する。

本実施例では、図１３のように該界面と各格子との切片をサンプリング点とすること各線分の中心点もサンプリング点として採用することにしている。図１３（a）が該物体の界面を表現する図形を表したものである。格子上に定義されている。図１３（b）において、該図形と格子との切片（白抜き丸）を計算する。この切片で、図１３（ｃ）のように該図形を近似する。（本実施例では、この新たな図形を与えられた界面であると認識している。）この近似された図形に対して、図１３（ｄ）のようにその中心点をもサンプリング点としている。つまり、図１３（ｄ）の白抜きの丸が、本実施例のサンプリング点である。本実施例ではこのようにしたが、もちろん図１３（ｃ）の近似を行わずに元の図形自身上の点を中間のサンプル点としてもよい。更には、格子スケールで適当に元の図形の表す界面自身に格子間隔と同じ程度となるように点を割り振ってもよい。例えば、擬似乱数を利用してもよい。

これらのサンプル点をステップＳ５では、点ｐからの距離をs[i] (i=1, ‥‥ Ｍ)とすると、d(p,Γ)の格子間隔の約２０分の１程度の近似精度をもつ近似距離Ｄ(p, Γ)を
Ｄ(p, Γ) ：= min d (p,s[i])
全てのiについて
として定義する。

また上述の通り、近似距離方向ベクトルは、格子点を始点、始点における近似距離の定義に採用した始点とペアを成すサンプリング点を終点とするベクトルとして定義する。

ここで、格子間隔の約２０分の１の精度を持つというのは1)通常、格子上で場を定義し、それらの離散的値によって物理現象を模擬する場合、距離を求めることが必要となる点が格子点あるいは双格子空間の格子点上に限られることと２）図１０の図形の場合を考えれば理解できる。図１０は格子点上で物体の界面からの距離を定義することを想定し、上記アルゴリズムによる物体の界面上のサンプリング点からの距離と、物体の界面からの距離の誤差が最悪の場合を考察している。図１０において、Δｘ、Δｙで定義した正方格子を考えている。Ｌ（ｘ）が物体の界面から格子点への距離を表している。他方、L1, L0がそれぞれ、近似距離である。界面の格子のｘ方向の切片をｘとして、ｘをパラメータとして本発明の近似距離であるmin(L0, L1)が真の距離Ｌ(x)とどのくらい異なるかを示したものである。誤差εとして、
ε = （min(L0, L1) - L(x)）/L(x)
を採用した。図１０（a）がその際の値を％表示したものである。ｘをパラメータとして横軸を動かした際に、最大の誤差が４%以下になっていることから、格子点のみで距離を計算する場合の最悪の近似が２０分の１程度で抑えられることを示している。

これにより、平方根という代数的関数を各サンプリング点に応じて１回のコールとそれらの最小を求める操作で本発明で利用する近似距離が求まる。本実施例では２次元の場合を取り扱ったが３次元の場合も同様であり、３次元の場合は浮動小数点演算の最大保証値までを利用した計算精度で求めようとすると膨大な計算量が必要となるのに対して、本アルゴリズムを採用すると１０倍程度の計算量の低減が実現できる。

本実施例によれば、これらのサンプリング点を格子との切片を利用して作成するために、格子間隔によってその精度によってスケーリングされている。

格子上で場を定義し、それらの離散的値によって物理現象を模擬する場合、重要となるのはこの格子間隔という計算における精度保証である。格子間隔を小さくすればする程、漸近的に、実際の現象との差（例えば、sup-normをとっても良い）が小さくなることが期待される。数値計算において、この格子間隔を単位とした漸近性が重要である。が、従来の計算方法では、これらの精度保証せず、計算機の浮動小数点の精度分の距離を計算していた。格子上で場を定義し、それらの離散的値によって物理現象を模擬する事を目的とした場合、これは明らかに不必要な精度（オーバースペック）であり、そのために計算速度を低下させていることは改善されなければならない。

本発明はこのような観点から、上記に示したように、距離計算に精度保証という概念を幾何学的に持ち込み、これを接触判定に応用することにより計算速度の向上を目指したものである。

また、本実施例では格子の格子点を代表点として物体の距離を計算して接触判定に用いたが、胞中心を代表点とした場合も同様に求まる。胞中心について直接、近似距離を計算してもよいし、各格子点での距離を利用して胞内の点の距離を補間したものを接触判定に用いてもよい。

また、本実施例の格子との切片を利用したサンプリング方法は、提案する近似距離に格子間隔という近似単位が持ち込まれる。このことにより、例えば図１４に示すような適合格子や非構造格子を用いた場合に関しても、細かい部分は小さな近似精度で、大まかでよい部分は大まか精度で距離を近似できる。よって、本発明の方法は正方構造格子のみではなく、その他の格子系においても適用でき、十分な効果を得ることができる。

（実施例２）
本発明の符号付き距離関数の計算方法を３次元の図形に対して適用した例について述べる。

空間に３次元の物体Ａ，Ｂを含む対象空間があるとする。続いて、対象空間を所望の空間解像度で胞を形成するように領域分割をし、空間固定格子を生成する。２次元の場合の胞は正方形、あるいは長方形などの多角形であったが、実施例１の２次元の場合と違って、３次元では胞は直方体、立方体などの多面体になる。ここでは、規則的に格子を形成したとして、胞は直方体、もしくは、立方体であるとする。さらに、物体Ａ，Ｂの界面が横切る各胞において、胞を構成する辺と物体との交点を胞内で結んだ面分を定義し、その面分をつなぎ合わせて物体形状を多面体図形（以降、近似多面体と呼ぶ。）により近似する。次に各格子点において、２次元の場合と同様にして符号付き近似距離関数を定義する。このとき、各格子点と近似多面体の各面との距離を計算するが、実施例１と同様の観点により、符号付き距離関数の値として、実際の距離を使用するのではなく、近似多面体と格子との交点をサンプリング点として、格子点からサンプリング点のまでの距離を近似距離として使用する。

さらに実施例１と同様にして、ペナルティ法により接触応力を計算する。本実施例においては接触応力は面応力となるため、近似多面体を構成する各面に対して接触応力を面積積分して接触力を求める。

実施例１の２次元場合と比べ、幾何学的な構成要素は変わるが、一連の処理の考え方は２次元の場合とほぼ同じであるため、本実施例を実施するための処理は実施例１の処理手順にわずかな修正を加えるだけで実現できる。

また、本実施例では実施例１と同様に、格子の格子点を代表点として物体の距離を計算して接触判定に用いたが、胞中心を代表点とした場合も同様に求まる。胞中心について直接、近似距離を計算してもよいし、各格子点での距離を利用して胞内の点の距離を補間したものを接触判定に用いてもよい。

また、実施例１の場合と同様に適合格子や非構造格子を用いた場合に関しても、細かい部分は小さな近似精度で、大まかでよい部分は大まか精度で距離を近似できることは自明である。よって、本発明の方法は正方構造格子のみではなく、その他の格子系においても適用でき、十分な効果を得ることができる。

本発明を計算処理装置を用いて実施する際の手順を示すフローチャート。実施例１の２次元空間内の物体Ａ，Ｂを示した図。実施例１において処理の手順を示す図。実施例１において処理の手順を示す図。実施例１において処理の手順を示す図。実施例１において処理の手順を示す図。実施例１において処理の手順を示す図。実施例１において処理の手順を示す図。本発明において線分との距離の計算方法の説明図。近似距離の精度に関する説明図。本発明の接触判定方式を搭載した処理装置の説明図。実施例１において処理の手順を示す図。本実施例１のサンプリング点の決定方法の説明図。本発明の適合格子への適用に関する説明図。実施例１において処理の手順を示す図。実施例１において処理の手順を示す図。実施例１において処理の手順を示す図。

符号の説明

１１１計算機の本体
１１２ディスプレイ
１１３フレキシブルディスクドライブ
１１４マウス
１１５キーボード
１１０フレキシブルディスク
１１６ＣＰＵ
１１７ハードディスク
１１８メモリ
Ｅ直交固定格子
Ｓａ物体Ａの形状を格子Ｌａの解像度で近似した折れ線図形
Ｓｂ物体Ｂの形状を格子Ｌｂの解像度で近似した折れ線図形
Ｓａ' 物体Ａの形状を格子Ｅの解像度で近似した折れ線図形
Ｓｂ' 物体Ｂの形状を格子Ｅの解像度で近似した折れ線図形
ｆａ物体Ａの符号付き近似距離関数
ｆｂ物体Ａの符号付き近似距離関数

Claims

二次元または三次元内の格子上に定義された第一の物体、第二の物体の物理量を数値的に模擬する際に必要となる接触判定、及び接触力を算出する方法において、格子間隔でパラメータ化された近似距離を計算する近似距離計算手段を有し、第一の物体、第二の物体それぞれの表面上に取った有限個の代表点同士を第一の物体、第二の物体それぞれについて幾何学要素でつなぎ合わせた近似図形により第一の物体、第二の物体それぞれの表面形状を近似し、該第一の物体、第二の物体と該第一の物体、第二の物体が定義された格子の格子点との該近似距離計算手段を用いて計算した符号付き近似距離を該格子点において計算することにより所望の領域で符号付き近似距離関数を定義し、該近似図形を構成する各代表点における該符号付き近似距離関数の値により、接触判定、及び接触力算出を行うことを特徴とする接触判定方式。