(57)【特許請求の範囲】
【請求項１】
移動境界を有する流体現象を数値的に解析する流体シミュレーションで用いられ、固定格子上を移動する移動境界をレベルセット法を用いて設定する、移動境界の設定方法において、
移流方程式に従って移動する移動境界の設定に必要とされるレベルセット関数から独立した関数として、移動境界を表すゼロレベルセットを追跡することを目的とするゼロレベルセット追跡関数であって、そのゼロ等高面がゼロレベルセットであり、レベルセット関数と同一の形を有するゼロレベルセット追跡関数を準備するステップと、
前記ゼロレベルセット追跡関数を前記移流方程式に従って、当該移流方程式の計算を安定的に行うために要求される安定性の条件から決定される時間間隔分だけ時間発展させる時間発展ステップと、
前記時間発展ステップにより時間発展された後のゼロレベルセット追跡関数に基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数を構築する境界設定ステップとを含み、
前記時間発展ステップおよび前記境界設定ステップを、所定の終了時間に達するまで繰り返すことを特徴とする、移動境界の設定方法。
【請求項２】
前記境界設定ステップにより構築されたレベルセット関数に基づいて、移動境界の内外にある物質の識別に用いられる密度関数を求めるステップをさらに含むことを特徴とする、請求項１記載の方法。
【請求項３】
前記時間発展ステップにおいて、前記移流方程式をＣＩＰ法により解くことを特徴とする、請求項１または２記載の方法。
【請求項４】
前記固定格子は直交固定格子であることを特徴とする、請求項１乃至３のいずれか記載の方法。
【請求項５】
移動境界を有する流体現象を数値的に解析する流体シミュレーションで用いられ、固定格子上を移動する移動境界をレベルセット法を用いて設定する、移動境界の設定方法を実現するコンピュータプログラムにおいて、
コンピュータに対して実行させる手順として、
移流方程式に従って移動する移動境界の設定に必要とされるレベルセット関数から独立した関数として、移動境界を表すゼロレベルセットを追跡することを目的とするゼロレベルセット追跡関数であって、そのゼロ等高面がゼロレベルセットであり、レベルセット関数と同一の形を有するゼロレベルセット追跡関数を準備する手順と、
前記ゼロレベルセット追跡関数を前記移流方程式に従って、当該移流方程式の計算を安定的に行うために要求される安定性の条件から決定される時間間隔分だけ時間発展させる時間発展手順と、
前記時間発展手順により時間発展された後のゼロレベルセット追跡関数に基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数を構築する境界設定手順とを含み、
前記時間発展手順および前記境界設定手順を、所定の終了時間に達するまで繰り返しコンピュータに対して実行させることを特徴とする、移動境界の設定方法を実現するコンピュータプログラム。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、移動境界を有する流体現象等を数値的に解析する手法に係り、とりわけ、固定格子上を移動する移動境界をレベルセット法を用いて設定する、移動境界の設定方法およびそれを実現するコンピュータプログラムに関する。
【０００２】
【従来の技術】
複雑な移動境界（気液、液液および気液固等）を有する液体現象は、我々の身の回り（例えば、河川やキッチン、体内等）の至るところに現れる。逆に、実際の流体現象としては、複雑な移動境界を有しない流体現象の方がむしろ少ないと言える。
【０００３】
このような流体現象に対する研究は、物理学や工学、医学等の様々な分野で重要であると考えられており、実験的な手法による解析に加えて、数値的な手法による解析も求められている。
【０００４】
例えば河川の浸食や堤防の破損の原因となる激流を減衰させるための構造物（ダム等）の設計は、大変形を伴う移動境界および複雑な移動境界の捕獲が考慮された適当な流体シュミレーション手法がないため大規模な実験や経験則のみによらざるを得ず、有効な手法の開発が望まれている。また、同様な課題は、海岸における突堤やテトラポットの設計および配置、工業装置における冷却パイプの振動破断を防止するためのパイプ形状の設計等にも共通するものである。また、有効な手法が得られればインクジェットやコーティング、洗浄等の効率化にも応用が可能である。
【０００５】
しかしながら、複雑な移動境界を有する流体現象を数値的な手法により解析する場合には、境界条件の扱いが複雑であり、また現象のスケールと複雑な移動境界のスケールとを同時に扱う必要があること等から、現象の解析には様々な困難を伴う。
【０００６】
複雑な移動境界を有する流体問題に対しては、固定格子（Euler grid）を用いることが有効である。これは、固定格子によれば、複雑な境界形状であっても容易に表現することができ、また境界の移動に伴ってその都度格子を生成する必要がないからである。また、固定されている格子点上を境界が移動するだけであるので、格子が捻れて計算が破綻するといったことがなく、境界の大変形やトポロジーの変化等にも容易に対応することができるからである。
【０００７】
しかしながら、固定格子を用いる場合には、境界での数値拡散やメッシュの影響等の問題が当然に起こり得る。そのため、様々な工夫が必要となる。
【０００８】
現在、このような問題に対処するための手法として、レベルセット法（level set method）が有効であると考えられている（文献１（S. Osher and J.A. Sethian, J. Comput. Phys. 79, 12 (1988)）および文献２（J.A. Sethian, Level Set Methods and Fast Marching Methods, Cambridge University Press (1999)）参照）。なお、レベルセット法は、現在までに様々な問題に適用されて成果を挙げている（上記文献２および文献３（K. Yokoi and F. Xiao, Phys. Rev. E 61 R1016 (2000)）参照）。
【０００９】
レベルセット法は、移動境界を設定（捕獲）するための手法の一つであり、（Ｎ−１）次元の境界をＮ次元の関数のゼロ等高面（zero-contour）（値が０である等値線や等値面）を用いて表現する。このＮ次元の関数をレベルセット関数（level set function）と呼び、そのゼロ等高面をゼロレベルセット（zero-level set）と呼ぶ。この手法では、境界自体を直接追跡するのではなく、レベルセット関数を介して、そのゼロレベルセットとして求められる境界を間接的に追跡する。
【００１０】
ここで、レベルセット関数としては、次式（１）の条件を満たす関数ｆlsが用いられる。
【００１１】
【数１】
0003978534.tif 000002
【００１２】
また、レベルセット関数ｆlsの時間発展は、次式（２）の移流方程式に従って計算される。なお、次式（２）において、ｕは流体の速度を表している。
【００１３】
【数２】
0003978534.tif 000003
【００１４】
ここで、図９に示すような１次元の流体問題を例にとると、上式（１）に示すレベルセット関数ｆlsは、符号２０に示すようなものとなり、液体と気体との境界２１，２２がレベルセット関数ｆlsのゼロ等高面（ゼロレベルセット）として求められる。
【００１５】
このため、上式（１）に示すレベルセット関数ｆlsを上式（２）の移流方程式に従って時間発展させ、ある時刻でのレベルセット関数ｆlsのゼロレベルセットを求めることにより、境界を追跡することができる。
【００１６】
なお、上式（１）および図９から分かるように、レベルセット関数ｆlsは傾きが１の関数であり、ある点２３におけるｆlsの値ａは、点２３の境界２１からの距離ｂと等しい。すなわち、レベルセット関数ｆlsは、境界からの距離を表す距離関数としての性質を有している。
【００１７】
このようなレベルセット法では、境界を直接追跡することなく、レベルセット関数ｆlsを用いて間接的に追跡するので、境界のトポロジーが変化するような場合（例えば、図１０(a)(b)(c)に示すように、２つの水滴（液体）が衝突するような場合）でも、特別な処理を行うことなく自然に対応することができる。また、直交固定格子上で曲面からなる境界を精度良く捕らえるには、境界をある程度スムースにする必要があるが、レベルセット法では、境界のスムーシングの幅を容易にコントロールすることができるので、曲面を精度良く捕らえることができる。さらに、表面張力の計算等に必要な曲率κの計算には、境界に対する単位法線ベクトルの発散を求める必要があり、境界の周辺での単位法線ベクトルが必要になるが、レベルセット法では、レベルセット関数ｆlsを用いることにより、次式（３）に従って曲率κを精度良く計算することができる。
【００１８】
【数３】
0003978534.tif 000004
【００１９】
なお、このようにして求められた曲率κに基づいて、ＣＳＦ（Continum Surface Force）モデル（文献４（J.U. Brackbill, D.B. Kothe and C. Zemach, J. Comput. Phys. 100, 335 (1992)）参照）等を用いて、容易に表面張力を計算することができる。
【００２０】
【発明が解決しようとする課題】
上述したように、レベルセット法は、移動境界を設定（捕獲）する手法として優れたものであるが、移動境界の設定に必要とされるレベルセット関数ｆlsには一般に、上式（２）の移流方程式に従って時間発展させるにつれてその関数の性質（特に、距離関数としての性質）が損なわれていくという問題がある。
【００２１】
すなわち、ある時刻でのレベルセット関数ｆlsが図１１(a)に示すように距離関数の性質を満たしていても、上式（２）の移流方程式により時間発展させた後のレベルセット関数ｆlsは、図１１(b)に示すように一般に距離関数の性質を満たさなくなる。
【００２２】
このため、従来のレベルセット法では、時間発展させた後のゼロレベルセットが境界を正しく表現しているものと仮定し、このゼロレベルセットを基準としてその他の部分に対して距離関数の性質を回復させるための処理、すなわち再初期化処理（re-initialization）を行っている。
【００２３】
しかしながら、このような再初期化処理に伴って、ゼロレベルセットも若干ずれることとなるので、これにより境界の位置に関してエラーが混入する（図１１(c)）。このエラーは１計算ステップ分で見れば僅かであるが、このエラーが混入したレベルセット関数（図１１(d)）を用いて次回以降の計算ステップを繰り返すこととなるので、このようなエラーが毎計算ステップごとに蓄積してしまう。このため、時間発展させるにつれて境界の位置が正しい位置から大きくずれてしまう。このようなエラーの蓄積は、流体問題の種類によっては深刻な問題となり、特に複雑な移動境界を有する流体問題の場合には、このエラーは無視できないものとなる。
【００２４】
本発明はこのような点を考慮してなされたものであり、レベルセット法において問題となる再初期化処理時のエラーの蓄積を効果的に防止して、固定格子上で移動境界を精度良く捕らえることができる、移動境界の設定方法およびそれを実現するコンピュータプログラムを提供することを目的とする。
【００２５】
【課題を解決するための手段】
本発明は、その第１の解決手段として、移動境界を有する流体現象を数値的に解析する流体シミュレーションで用いられ、固定格子上を移動する移動境界をレベルセット法を用いて設定する、移動境界の設定方法において、移流方程式に従って移動する移動境界の設定に必要とされるレベルセット関数から独立した関数として、移動境界を表すゼロレベルセットを追跡することを目的とするゼロレベルセット追跡関数であって、そのゼロ等高面がゼロレベルセットであり、レベルセット関数と同一の形を有するゼロレベルセット追跡関数を準備するステップと、前記ゼロレベルセット追跡関数を前記移流方程式に従って、当該移流方程式の計算を安定的に行うために要求される安定性の条件から決定される時間間隔分だけ時間発展させる時間発展ステップと、前記時間発展ステップにより時間発展された後のゼロレベルセット追跡関数に基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数を構築する境界設定ステップとを含み、前記時間発展ステップおよび前記境界設定ステップを、所定の終了時間に達するまで繰り返すことを特徴とする、移動境界の設定方法を提供する。
【００２６】
なお、上述した第１の解決手段においては、前記境界設定ステップにより構築されたレベルセット関数に基づいて、移動境界の内外にある物質の識別に用いられる密度関数を求めるステップをさらに含むことが好ましい。また、前記時間発展ステップにおいて、前記移流方程式をＣＩＰ法により解くことが好ましい。さらに、前記固定格子は直交固定格子であることが好ましい。
【００２７】
本発明は、その第２の解決手段として、移動境界を有する流体現象を数値的に解析する流体シミュレーションで用いられ、固定格子上を移動する移動境界をレベルセット法を用いて設定する、移動境界の設定方法を実現するコンピュータプログラムにおいて、コンピュータに対して実行させる手順として、移流方程式に従って移動する移動境界の設定に必要とされるレベルセット関数から独立した関数として、移動境界を表すゼロレベルセットを追跡することを目的とするゼロレベルセット追跡関数であって、そのゼロ等高面がゼロレベルセットであり、レベルセット関数と同一の形を有するゼロレベルセット追跡関数を準備する手順と、前記ゼロレベルセット追跡関数を前記移流方程式に従って、当該移流方程式の計算を安定的に行うために要求される安定性の条件から決定される時間間隔分だけ時間発展させる時間発展手順と、前記時間発展手順により時間発展された後のゼロレベルセット追跡関数に基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数を構築する境界設定手順とを含み、前記時間発展手順および前記境界設定手順を、所定の終了時間に達するまで繰り返しコンピュータに対して実行させることを特徴とする、移動境界の設定方法を実現するコンピュータプログラムを提供する。
【００２８】
本発明によれば、レベルセット法で用いられるレベルセット関数ｆlsから独立した関数としてゼロレベルセット追跡関数を準備し、このゼロレベルセット追跡関数を移流方程式に従って時間発展させ、ゼロレベルセット追跡関数に基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数を構築しているので、時間発展の計算によりレベルセット関数にエラーが蓄積することがなく、従来のレベルセット法で問題となっていた再初期化処理時のエラーの蓄積の問題を効果的に防止することができる。
【００２９】
【発明の実施の形態】
以下、図面を参照して本発明の実施の形態について説明する。ここで、図１は本発明による移動境界の設定方法の一実施の形態を説明するためのフローチャート、図２(a)(b)(c)(d)は本実施の形態の概要を１次元の移動境界問題を例に挙げて説明するための概念図である。なお、本実施の形態に係る移動境界の設定方法は、直交固定格子上を移動する移動境界をレベルセット法を用いて設定するものであり、移動境界を有する流体現象等の解析に用いられる。なお、固定格子としては、直交固定格子以外にも、三角格子等の任意の固定格子を用いることができる。
【００３０】
図１に示すように、本実施の形態に係る移動境界の設定方法においては、まず、レベルセット法で用いられるレベルセット関数ｆlsから独立した関数として、移動境界を表すゼロレベルセットを追跡することのみを目的とするゼロレベルセット追跡関数ｆzeroを準備する（ステップ１０１）。なお、広い意味でいえば、ゼロレベルセット追跡関数ｆzeroもレベルセット関数であり、そのゼロ等高面はゼロレベルセットであるが、本明細書では、距離関数の性質を持つ関数のみをレベルセット関数と呼び、このようなレベルセット関数のゼロ等高面のみをゼロレベルセットと呼ぶことにする。
【００３１】
次に、ゼロレベルセット追跡関数ｆzeroを次式（４）の移流方程式に従って、安定性の条件から決定される時間間隔分だけ時間発展させる（ステップ１０２）。具体的には例えば、図２(a)(b)に示すように、ゼロレベルセット追跡関数ｆzeroを境界の移動に合わせて右方向に移動させる。なおこのとき、ゼロレベルセット追跡関数ｆzeroは時間発展に伴ってその形状が崩れる。
【００３２】
【数４】
0003978534.tif 000005
【００３３】
なお、上式（４）の移流方程式の解法には、高精度でかつロバストな方法であるＣＩＰ法（Cubic Interpolated Propagation Method）（文献５（T. Yabe and T. Aoki, Comput. Phys. Commun. 66, 219 (1991)）および文献６（T. Nakamura and T. Yabe, Comput. Phys. Commun. 120, 122 (1999)）参照）を用いることが好ましく、これにより、ゼロレベルセット追跡関数ｆzeroの形状の崩れを抑えることができる。なお、移流方程式の解法には、ＣＩＰ法以外にも、ＰＰＭ（Piecewise Parabolic Method）法等の任意の方法を用いることができる。
【００３４】
その後、時間発展された後のゼロレベルセット追跡関数ｆzeroに基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数ｆlsを構築する（ステップ１０３）。具体的には例えば、図２(c)に示すように、ゼロレベルセット追跡関数ｆzeroのゼロ等高面に基づいて、移動境界を表すゼロレベルセットを求め、この求められたゼロレベルセットに基づいて、距離関数の性質を持つレベルセット関数ｆlsを構築する。なお、レベルセット関数ｆlsは、境界のスムーシングや曲率の計算等のために用いられる。
【００３５】
以下、ステップ１０３における処理（ゼロレベルセット追跡関数ｆzeroからレベルセット関数ｆlsを構築する処理）の詳細について説明する。
【００３６】
ステップ１０３においては、まず、ゼロレベルセット追跡関数ｆzeroのゼロ等高面を、各格子点間を補間することにより求める。
【００３７】
次に、このようにして求められたゼロ等高面を、移動境界を表すゼロレベルセットとし、境界から１セル以内の格子点のレベルセット関数ｆlsを、ファーストマーチング法（fast marching method）（上記文献２および文献７（D. Adalsteinsson and J.A. Sethian, J. Comput. Phys. 148, 2 (1999)）参照）に従って、次式（５）のEikonal方程式を計算することにより求める。
【００３８】
【数５】
0003978534.tif 000006
【００３９】
次に、このようにして求められたレベルセット関数ｆlsの値を固定し、Sussmann等の反復計算による方法（文献８（M. Sussman, P. Smereka and S. Osher, J. Comput. Phys. 114, 146 (1994)）参照）に従って、次式（６）が収束するまで収束計算を行うことにより、１セル以遠の部分のレベルセット関数ｆlsを構築する。なお、Sussman等の方法では一般に、反復計算の間にゼロレベルセットがずれるという問題があるが、ここでは、ゼロレベルセットの周辺の値がファーストマーチング法により事前に計算され、実際の収束計算時には固定されているため、その問題は起こらない。
【００４０】
【数６】
0003978534.tif 000007
【００４１】
ここで、τは反復計算のための時間であり、実時間と同じである必要はない。また、Ｓは平滑化符号関数（smoothed sign function）であり、次式（７）により表される。
【００４２】
【数７】
0003978534.tif 000008
【００４３】
なお、ゼロレベルセットの周辺のレベルセット関数ｆlsを固定しているため、計算の安定性が問題とならない場合には、Ｓ（ｆls）＝１とすることも可能である。
【００４４】
ここで、上式（６）による収束計算における初期値としては、前の計算ステップのｆlsの値や、１、−１等を用いるのではなく、上式（２）の移流方程式により計算した値を用いるようにしてもよく、これにより、収束に要する時間を大幅に短縮することができる。このような収束計算には、上述したＣＩＰ法を用いることができる他、１次の風上差分法を用いることができる。これは、１次の風上差分法のような比較的低精度の方法を用いても、その方法により混入されたエラーが上式（６）による収束計算により補正されるからである。なお、計算時間に拘らない場合には、上式（６）による収束計算ではなく、上述したファーストマーチング法を用いることもできる。
【００４５】
その後、このようにして構築されたレベルセット関数ｆlsに基づいて、移動境界の内外にある物質の識別に用いられる密度関数（color (density) function）ｆcolorを求める（ステップ１０４）。なお、密度関数ｆcolorは、次式（８）（９）に従って、レベルセット関数ｆlsから平滑化ヘビサイド関数（smoothed Heaviside function）Ｈαを用いて計算することができる。ここで、２αは境界のスムーシングの幅である。
【００４６】
【数８】
0003978534.tif 000009
【００４７】
そして、所定の終了時間に達するまでステップ１０２〜１０４における処理を繰り返す（ステップ１０５）。なお、ステップ１０２に戻って処理を続ける場合、時間発展の計算の基礎となる関数はゼロレベルセット追跡関数ｆzeroである。すなわち、ステップ１０３において、ゼロレベルセット追跡関数ｆzeroに基づいてレベルセット関数ｆlsを構築するときには、従来のレベルセット法と同様にエラーが混入することとなるが（図２(c)参照）、ステップ１０２における時間発展の計算には、エラーが混入していないゼロレベルセット追跡関数ｆzeroが用いられる（図２(d)参照）。
【００４８】
このように本実施の形態によれば、レベルセット法で用いられるレベルセット関数ｆlsから独立した関数としてゼロレベルセット追跡関数ｆzeroを準備し、このゼロレベルセット追跡関数ｆzeroを移流方程式に従って時間発展させ、ゼロレベルセット追跡関数ｆzeroに基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数ｆlsを構築しているので、時間発展の計算によりレベルセット関数ｆlsにエラーが蓄積することがなく、従来のレベルセット法で問題となっていた再初期化処理時のエラーの蓄積の問題を効果的に防止することができる。
【００４９】
なお、上述した実施の形態に係る移動境界の設定方法は、図３に示すようなコンピュータシステム１０上で稼働するコンピュータプログラムとして実現することができる。ここで、コンピュータシステム１０は、バス１１と、バス１１に接続されたプロセッサ１２、メモリ１３およびハードディスク１４と、バス１１に接続された周辺機器（キーボードやマウス等の入力装置１５、ディスプレイやプリンタ等の出力装置１６、ＦＤドライブ１７およびＣＤ−ＲＯＭドライブ１８）とを備えている。そして、上述したようなコンピュータプログラムは、フレキシブルディスク１９ａおよびＣＤ−ＲＯＭ１９ｂ等のようなコンピュータ読み取り可能な記録媒体に格納され、メモリ１３やハードディスク１４等に読み込まれる。そして、コンピュータシステム１０において、プロセッサ１２から逐次読み出されて実行されることにより上述したような手順を実現することができる。
【００５０】
【実施例】
次に、上述した実施の形態の具体的実施例について述べる。
【００５１】
まず、上述した実施の形態をZalesak問題（文献９（S.T. Zalesak, J. Comput. Phys. 31, 335 (1979)）参照）に適用した場合の実施例について述べる。
【００５２】
図４に示すように、本実施例の移動境界問題は、一部に切り欠きを有する円形状の二次元の剛体３１を矢印３２に沿って反時計方向に回転させるものである。
【００５３】
ここで、剛体３１の初期条件は、次式（１０）とした。
【００５４】
【数９】
0003978534.tif 000010
【００５５】
また、計算のための格子は１００×１００の直交固定格子とし、剛体３１の回転速度は、次式（１１）のような速度場の下で、１回転の計算を６２８回の計算ステップで行うようなものとした。（すなわち、時間間隔Δｔを１／６２８とした。）
【００５６】
【数１０】
0003978534.tif 000011
【００５７】
さらに、レベルセット関数ｆlsは図５(a)に示すものとした。なお、このようなレベルセット関数ｆlsに対応する、初期状態でのゼロレベルセット（境界）および密度関数ｆcolorはそれぞれ、図５(b)(c)に示すようなものとなる。
【００５８】
このような条件の下で、レベルセット関数ｆlsと同一の形を有するゼロレベルセット追跡関数ｆzeroを準備し、このゼロレベルセット追跡関数ｆzeroを上式（１１）に従ってＣＩＰ法により時間発展させ、ゼロレベルセット追跡関数ｆzeroに基づいて、移動境界を表すゼロレベルセットを求めるとともに、この求められたゼロレベルセットに基づいて、移動境界の設定に必要とされるレベルセット関数ｆlsを構築し、移動境界を設定した。
【００５９】
図６(a)(b)および図６(c)(d)はそれぞれ、剛体３１の１回転後および１０回転後の境界の様子を示す図である。なお、図６(a)(c)はレベルセット関数ｆlsに基づいて求められたゼロレベルセット（境界）を厳密解（初期条件）とともに示す図であり、図６(b)(d)はレベルセット関数ｆlsに基づいて求められた密度関数ｆcolorを示す図である。
【００６０】
図６(a)(b)(c)(d)に示すように、本実施例では、剛体３１の１回転後はもちろん１０回転後であっても境界を精度良く捕らえていることが分かる。
【００６１】
ここで、比較例として、上述した従来のレベルセット法を用いて剛体３１を１回転させた場合の境界の様子を図７(a)(b)に示す。図７(a)はゼロレベルセット（境界）を厳密解（初期条件）とともに示す図であり、図７(b)は密度関数ｆcolorを示す図である。なお、この比較例におけるレベルセット法では、レベルセット関数ｆlsを時間発展させる際の移流方程式の解法としてＣＩＰ法を用いた。
【００６２】
図７(a)(b)に示すように、この比較例では、１回転後にもかかわらず、移動境界が初期状態からかなりずれていることが分かる。なお、このずれは時間とともに拡大した。
【００６３】
次に、上述した実施の形態を流体問題に適用した場合の実施例について述べる。
【００６４】
図８(a)(b)に示すように、本実施例の流体問題は、長方形の剛体３３を水面３４に落下させるものである。なお、計算のための格子は４０×４０の直交固定格子とし、時間間隔（Δｔ）は安定性の条件を満たすように決定した。また、流体シミュレーションの手法としては、圧力のポアソン方程式をベースにした混相流の計算手法であるＣＵＰ法（文献１０（T. Yabe and P.Y. Wang, J. Phys. Soc. Jpn. 60, 2105 (1991)）参照）を用いた。
【００６５】
その結果、図８(b)に示すように、本実施例においては、比較的少ない格子数にもかかわらず、剛体３３の形状を良く捕らえていることが分かる。これに対し、図８(c)に示す比較例（従来のレベルセット法を用いたもの）では、再初期化処理時のエラーにより、剛体３３の角が丸くなる傾向があることが分かる。
【００６６】
【発明の効果】
以上説明したように本発明によれば、レベルセット法において問題となる再初期化処理時のエラーの蓄積を効果的に防止して、固定格子上で移動境界を精度良く捕らえることができる。従って、本発明は、河川の激流や、海岸の大規模波動、工業装置における流体現象等に対応する構造物の設計等に極めて有用である。
【図面の簡単な説明】
【図１】本発明による移動境界の設定方法の一実施の形態を説明するためのフローチャート。
【図２】図１に示す実施の形態の概要を１次元の移動境界問題を例に挙げて説明するための概念図。
【図３】図１に示す実施の形態を実現するコンピュータプログラムが実行されるコンピュータシステムの一例を示す図。
【図４】本発明の一実施例としてのZalesak問題（二次元の剛体の回転問題）の概要を説明するための模式図。
【図５】図４に示す実施例におけるレベルセット関数、ゼロレベルセットおよび密度関数の一例を示す図。
【図６】図４に示す実施例における剛体の回転後（１回転後および１０回転後）の境界の様子を示す図。
【図７】従来のレベルセット法を用いた場合の剛体の回転後（１回転後）の境界の様子を示す図。
【図８】本発明の一実施例としての流体問題（長方形の剛体を水面に落下させる問題）の計算結果を示す図。
【図９】１次元の流体問題で用いられるレベルセット関数の一例を示す図。
【図１０】境界のトポロジーが変化するような１次元の流体問題の一例を示す図。
【図１１】従来のレベルセット法の概要を説明するための概念図。
【符号の説明】
ｆzero  ゼロレベルセット追跡関数
ｆls  レベルセット関数
１０  コンピュータシステム
２０  レベルセット関数
２１，２２  境界
【図１】
0003978534.tif 000012
【図２】
0003978534.tif 000013
【図３】
0003978534.tif 000014
【図４】
0003978534.tif 000015
【図５】
0003978534.tif 000016
【図６】
0003978534.tif 000017
【図７】
0003978534.tif 000018
【図８】
0003978534.tif 000019
【図９】
0003978534.tif 000020
【図１０】
0003978534.tif 000021
【図１１】
0003978534.tif 000022