JP3566009B2

JP3566009B2 - Particle type fluid simulation method and apparatus

Info

Publication number: JP3566009B2
Application number: JP34936796A
Authority: JP
Inventors: 滋大木; 直子松田; 裕之中; 達男横塚
Original assignee: Panasonic Corp; Matsushita Electric Industrial Co Ltd
Current assignee: Panasonic Corp; Panasonic Holdings Corp
Priority date: 1996-12-27
Filing date: 1996-12-27
Publication date: 2004-09-15
Anticipated expiration: 2016-12-27
Also published as: JPH10185755A

Description

【０００１】
【発明の属する技術分野】
本発明は、流体（液体又は気体）の形状の変化、又は流体から固体までの表面形状（界面状態）を定量的に表現し、動的及び静的な流体挙動を予測する粒子型流体シミュレーション方法及び装置に関するものである。
【０００２】
【従来の技術】
従来の流体解析技術は、２種類の手法に分類でき、１つは、流体の基本方程式であるナビエ・ストークスの方程式（Ｎａｖｉｅｒ−Ｓｔｏｋｅｓの式）を解くものと、もう１つは、自由表面を扱うＶＯＦ法（ＡＶｏｌｕｍｅｏｆＦｌｕｉｄＡｌｇｏｒｉｔｈｍ）である。前者は、流体の表面形状には触れておらず、流体内部の流速や圧力や流体抵抗について、計算するものであり、後者は、流体の表面形状に関して計算可能だが、正確な計算ができていないといった問題点が存在する。
【０００３】
具体的な説明を図１９〜図２１で行う。まず、前者の手法では、液体そのものが経時的に変化しても、そのような変化については計算していない。図１９（Ａ）から（Ｂ）に示すように液体６２を入れた容器６０を斜めに傾けた際に容器６０の傾斜に応じて液体全体の形状が変化するにもかかわらず、この種の問題については、取り扱えない問題が存在する。図中、６３は容器６０内の液体６２の界面である。ただし、図２０のように、液体が管全体の中を流れる速さｖや流跡や圧力差ΔＰや流体抵抗を解く計算は行えるが、流体の表面について形状を予測したり、液体の分離及び合体のような表面形状そのものに関しては、取り扱うことができない。後者の手法についても、図２１を利用して説明すると、図２１（Ｂ）のように分割された領域の中で体積の占有率で液体７２の形状を把握するように計算する方法では、容器７０内の液体７２の境界（液体７２の界面７１）については、図２１の（Ｃ）及び（Ｄ）に示すように分割された領域において体積占有率が同一であれば、どのような形状であるか正確に判断できない問題点が存在する。従って、液体７２の表面形状について正確に計算することができない。同様に、液体の分離及び合体のような表面形状の予測も取り扱うことができない。
【０００４】
しかし、近年、様々な液体を利用して、液体を固体表面上に正確かつ均質に膜を形成する技術が、半導体のみならず様々な技術分野に必要不可欠となっている。このとき、液体の表面形状がどのようになっているかを予測することは、非常に重要度を増してきている。例えば、図２２（Ａ）のように圧力Ｐの作用により容器９１から吐出される液体９２の全体の流速ｖのみならず、液ヌレやハジキに伴う容器９１の固体表面上と液体表面との関係で液体９２が横飛び９４を起こす現象は、従来の解析手法では、前述した理由で予測困難であるという問題点が存在する。図２２（Ｂ）において、９３は液ヌレにより液滴が容器９１の開口端面に付着した状態を示している。
【０００５】
【発明が解決しようとする課題】
従来の流体解析技術やＶＯＦ法に代表される自由表面解析技術においては、前述の具体例で説明したように、ミクロな液滴形状そのものの算出や液体と固体表面形状とのヌレ及びハジキに伴う液の横飛びした現象は、計算不可能であった。そのため、このような現象を再現するためには、試行錯誤の実験を繰り返す他なかった。
従って、本発明の目的は、従来の流体解析技術やＶＯＦ法に代表される自由表面解析技術が有していた上記の問題点を解決することにあって、流体のミクロな表面形状そのものの算出やヌレ及びハジキに伴う現象などをより正確に計算することが可能な粒子型流体シミュレーション方法及び装置を提供することである。
【０００６】
【課題を解決するための手段及びその作用効果】
上記目的を達成するために、本発明は以下のように構成する。
すなわち、本発明の第１態様にかかる粒子型流体シミュレーション方法によれば、流体を、連続体ではなく、各微小領域を仮想的な粒子とし、各粒子に働く力の相互作用を表現する粒子間の距離をパラメータとした関数に基づいて、各粒子の位置の経時的変化を決定し、流体の形状の変化を再現するものであって、
上記関数は、各粒子間の距離がある臨界点距離となるとお互いに斥力を生じ、この臨界点距離を越えて一定の距離になるとお互いに引力を生じることを特徴とする曲線又は直線からなり、この曲線又は直線は、各粒子の温度変化、速度変化、又は時間変化に伴い、曲線又は直線自体が変化することにより、非ニュートン性粘度の変化を表現し、流体の形状の変化を再現するようにしている。
本発明の第２態様にかかる発明では、第１態様において、前記関数に基づいて、各粒子の圧力の経時的変化をも決定し、流体の形状の変化を再現することもできる。
【０００７】
本発明の第３態様にかかる粒子型流体シミュレーション方法によれば、流体を、連続体ではなく、各微小領域を仮想的な粒子とし、各粒子と固体表面との関係において親水性及び撥水性を再現するために、各粒子と固体表面との力の相互作用を固体表面と粒子間の距離をパラメータとした関数に基づいて、各粒子の位置の経時的変化を決定し、流体と固体表面との親水性及び撥水性等の現象を再現するようにしている。
本発明の第４態様によれば、第３態様において、前記関数に基づいて、各粒子の圧力の経時的変化をも決定し、流体と固体表面との親水性及び撥水性等の現象を再現することもできる。
本発明の第５態様によれば、第３又は４態様において、前記関数は、各粒子と固体表面との間の距離がある臨界点距離となるとお互いに斥力を生じ、この臨界点距離を越えて一定の距離になるとお互いに引力を生じることを特徴とする曲線又は直線からなり、この曲線又は直線は、各粒子の温度変化、速度変化、又は時間変化に伴い、曲線又は直線自体が変化することにより、非ニュートン性粘度の変化を表現し、流体と固体表面との親水性及び撥水性等の現象を再現するようにすることもできる。
【００１０】
従って、本発明の上記種々の態様によれば、上記流体を仮想的な粒子に置き換えて表現し、各粒子間に働く力の相互作用を粒子間の距離をパラメータとした関数に基づいて、各粒子の位置、速度、密度、圧力、又は温度の経時的変化を決定することができる。
また、本発明の上記種々の態様によれば、上記粒子と上記固体表面との力の相互作用を固体表面と粒子間の距離をパラメータとした関数に基づいて、各粒子の位置、速度、密度、圧力、温度の経時的変化を決定することにより、上記流体と上記固体表面とのヌレ（親水性）及びハジキ現象（撥水性）を再現することができる。
また、各粒子間に働く力の相互作用を温度変化、速度変化、又は時間的変化に伴い各粒子間の距離をパラメータとした関数を変更制御することにより、非ニュートン性粘度のチクソ性等による時間的変化に伴い、粘度が変化する物性に対応することができる。
また、各粒子の質量、体積、温度、粘度、及び表面張力等のパラメータを複数持つことにより、２種類以上の異なる物性からなる流体、又は、２相流（固体及び液体、固体及び気体、液体及び気体、固体及び液体及び気体）以上の状態からなる物体の状態を再現することができる。
よって、本発明によれば、流体（気体又は液体）から固体までの表面形状を粒子型モデルで粒子間力を計算することで、流体の粘性や表面張力による濡れ性やミクロな表面形状などをより正確にシミュレーションすることができ、流体（液体又は気体）の形状の変化、又は流体から固体までの表面形状（界面状態）を定量的に表現し、動的及び静的な流体挙動を予測することができる。よって、例えば、液体塗布装置など流体吐出時の流体収納容器又は吐出ノズルなどの最適な形状や吐出流体の最適粘度を決定することができる。
【００１１】
【発明の実施の形態】
以下に、本発明にかかる実施形態を図面に基づいて詳細に説明する。
本出願は、本出願人が既に出願した特許出願平成８−１７８８９４号（発明の名称「液体及び濡れの評価装置」）は、流体（気体又は液体）を、連続体ではなく、各微小領域を仮想的な粒子として扱う手法について記述したものであり、具体的な計算方法に関してまでは詳細に記載したしたものではなかった。そこで、本発明者は、鋭意研究の結果、流体を粒子で表現した場合の具体的な計算のアルゴリズムを明確にすることができる発明を創作するに至ったものである。その内容は大略２つの点からなり、１点目は粒子間に働く力の相互作用を粒子間の距離をパラメータとした関数に基づき決定することと、２点目は、粒子と固体表面にかかる相互作用においても粒子と固体表面の最短距離をパラメータとした関数に基づき決定することで液体と固体表面との親水性及び撥水性の性質を再現する手法について述べるものである。
【００１２】
言い換えれば、本実施形態では、流体を、連続体ではなく、各微小領域を仮想的な粒子として扱い、これらの粒子のうちの各粒子に働く力の相互作用を粒子間の距離をパラメータとした関数に基づいて、各粒子の位置（これ以外に、例えば、速度、密度、圧力、温度など）の経時的変化を決定し、流体の形状の変化を再現するものである。そして、各粒子に働く力の相互作用を粒子間の距離をパラメータとした関数に基づいて、上記各粒子の圧力の経時的変化をも決定し、流体の形状の変化を再現することもできる。さらに、上記各粒子間に働く力の相互作用を表現する上記距離のパラメータの関数は、各粒子間の距離がある臨界点距離となるとお互いに斥力を生じ、この臨界点距離を越えて一定の距離になるとお互いに引力を生じることを特徴とする曲線又は直線からなり、この曲線又は直線は、各粒子の温度変化、速度変化、又は時間変化に伴い、曲線又は直線自体が変化することにより、非ニュートン性粘度のチクソ性による時間的変化に伴う粘度の変化を表現することもできる。
【００１３】
また、本実施形態では、さらに、流体を、連続体ではなく、各微小領域を仮想的な粒子として扱い、各粒子と固体表面との関係において親水性及び撥水性を再現するために、各粒子と固体表面との力の相互作用を固体表面と粒子間の距離をパラメータとした関数に基づいて、各粒子の位置（これ以外に、例えば、速度、密度、圧力、温度など）の経時的変化を決定し、流体と固体表面との親水性及び撥水性等の現象を再現するものである。そして、各粒子と固体表面との力の相互作用を固体表面と粒子間の距離をパラメータとした関数に基づいて、上記各粒子の圧力の経時的変化をも決定し、流体と固体表面との親水性及び撥水性等の現象を再現することもできる。さらに、上記各粒子と固体表面との間に働く力の相互作用を表現する上記距離のパラメータの関数は、各粒子と固体表面との間の距離がある臨界点距離となるとお互いに斥力を生じ、この臨界点距離を越えて一定の距離になるとお互いに引力を生じることを特徴とする曲線又は直線からなり、この曲線又は直線は、各粒子の温度変化、速度変化、又は時間変化に伴い、曲線又は直線自体が変化することにより、非ニュートン性粘度のチクソ性による時間的変化に伴う粘度の変化を表現することもできる。
また、上記各粒子には、各粒子の質量、体積、粘度、表面張力の粒子パラメータを複数個持つようにして、種々のモデルを表現することができるようにすることもできる。
【００１４】
以下、本発明の一実施形態にかかる粒子型流体シミュレーション方法及び装置について図面を参照しながら説明する。
図１は、上記本発明の上記実施形態にかかる粒子型流体シミュレーション方法を実施するためのフローチャートである。図２は、上記本発明の上記実施形態にかかる粒子型流体シミュレーション方法を実施するための上記シミュレーション装置のブロック図である。
図２において、キーボード又はマウスなどの入力部１から、解析すべき流体及び境界条件などの解析条件などを入力して記憶部２に記憶させたのち、演算部３で記憶部２に記憶されたデータを基に上記シミュレーション方法を実施し、シミュレーション結果を表示画面又は印刷装置などの出力部４で出力する。
【００１５】
上記演算部３で行う上記シミュレーション方法について、図１のフローチャートに基づいて説明する。
まず、ステップＳ１では、初期設定として流体の物性値、又は、流体及び固体表面の物性値を入力部１から記憶部２に設定入力して記憶させるとともに、初期配置として流体の位置及び形状、又は、流体及び固体表面の位置及び形状を同様に入力部１から入力して記憶部２で記憶させる。
次いで、ステップＳ２では、演算部３で、注目する粒子の近傍粒子を検索し、近傍粒子が存在する場合には上記注目粒子とその近傍粒子との中心間の距離計算を行うか、又は、注目する粒子の近傍固体表面を検索し、近傍固体表面が存在すれば注目粒子と近傍固体表面との距離計算を行う。
【００１６】
次いで、ステップＳ３では、上記演算部３での演算結果に基づき、入力部１で、粒子間又は粒子固体表面間の距離に依存したパラメータの選択を行う。粒子間同士の距離に依存したパラメータの選択を行う場合には、ステップＳ４で、液体の種類などに応じて粒子間同士の距離に依存した相互作用のパラメータの参照を行う。例えば、液体Ａと液体Ａとの間、液体Ｂと液体Ｂとの間、又は、液体Ａと液体Ｂとの間など解析すべき対象のそれぞれに応じて、温度変化に伴うパラメータ、速度変化に伴うパラメータ、時間変化に伴うパラメータなどのうちから適当なものを選択する。その後、ステップＳ６に進む。一方、上記演算部３での演算結果に基づき、ステップＳ３で粒子と固体表面との間の距離に依存したパラメータを選択する場合には、ステップＳ５に進む。このステップＳ５では、固体表面に応じて粒子と固体表面との間の距離に依存した相互作用のパラメータの参照を行う。例えば、固体表面Ｃ、又は、固体表面Ｄなど解析すべき対象のそれぞれに応じて、温度変化に伴うパラメータ、速度変化に伴うパラメータ、時間変化に伴うパラメータなどのうちから適当なものを選択する。その後、ステップＳ６に進む。なお、上記液体の種類などに応じて粒子間同士の距離に依存した相互作用の各種パラメータ、及び、上記粒子と上記固体表面との間の各種パラメータは記憶部２に予め記憶されており、適宜、入力部１からの入力操作により記憶部２のこれらのパラメータのデータを参照して選択する。
【００１７】
次いで、ステップＳ６では、上記選択されたパラメータを使用して演算部３で各粒子の受ける力の計算を行う。
次いで、ステップＳ６での計算結果に基づき、ステップＳ７において、演算部３により各粒子の位置の変更を行う。このとき、各粒子の位置以外に、又は、位置は変更することなく、圧力、速度、密度なども必要に応じて変化させる。
次いで、ステップＳ８で、演算部３により解析が終了か否かを判断し、終了でないならば、ステップＳ２に戻る。終了ならば、ステップＳ９で結果を出力部４により適宜出力して、解析を終了する。
【００１８】
ここで、上記ステップＳ６の各粒子の受ける力の計算及びステップＳ７の各粒子の位置の変更について、図３，４に基づいて、さらに、詳細に説明する。
図３，４は、上記実施形態にかかるステップＳ６，Ｓ７の動作において行われるアルゴリズムを記述したものである。ここで、各粒子は、大きく分けて以下の４つの力を受けて、現在の時刻から比べて、次の時刻において各粒子の位置、速度、圧力、温度、及び密度を変化させていると考える。
【００１９】
まず、１つ目の力は、粒子間同士に働く力で、この力により、流体内部の粘性や延性を表すことが可能となる。その計算方法は、粒子ｉと粒子ｊとを考える場合、粒子ｉと粒子ｊの中心間の粒子間距離ｒによって粒子同士間に生じる力Ｆが変化すると考える。これは、図５に示すように、一般に、２つの分子間では分子間力が作用することに鑑みて、粒子同士間においても粒子同士の距離をパラメータとした相互作用が働くものとして近似的に考えるものである。具体的には、縦軸を上記粒子同士間に生じる力Ｆとし、横軸を上記粒子間距離ｒ_ｘとした場合の両者の関係を図５にグラフで示す。このグラフが横軸を横切る点が以下に述べる臨海点距離ｒ_ｃである。よって、このグラフにおいて、粒子ｉと粒子ｊとの中心間の距離である粒子間距離ｒが、図５（Ａ）に示すように、ある臨界点距離ｒ_ｃより短くなると互いに斥力が生じ、図５（Ｂ）に示すように、粒子間距離ｒがある臨海点距離ｒ_ｃに等しくなると粒子間には力が働かなくなり、図５（Ｃ）に示すように、粒子間距離ｒがある臨界点距離ｒ_ｃより長くなると粒子間には互いに引き合う引力が生じる。さらに、図５には具体的には図示しないが粒子間距離がさらに大きくなると粒子間では互いに力を及ぼさなくなることがわかる。また、引力を表す面積が大きいほど基本的には粘度が大きくなり、非常に粘度が強くなると固体のように固まった物性になる。よって、この引力を表わす面積を変化させることにより、粘度を変化させることができる。また、臨界点の距離ｒ_ｃを変えることで、粒子の密度が変化することにもなる。すなわち、臨界点距離ｒ_ｃが大きくなれば粒子間の安定する距離が長くなるため粒子の密度が低くなる。これにより、粒子間の距離を変化させることにより、粒子の密度を変化させることが可能となる。
【００２０】
２つ目の力は、粒子と固体表面間に働く力で、この力により液体の親水性や撥水性を表現することが可能となる。その計算方法は、粒子間の力と同様に粒子の中心と固体表面との間の最短距離を計算し、その距離をパラメータとした関数として表現する。粒子同士の計算と基本的には、同様で、粒子と固体表面との相性によりその関数は変化する。
３つ目の力は、重力であり、全ての粒子に対して常に一定方向に重力が作用している。
４つ目の力は、その他の外力として、そのときにより、静電気力や慣性力など場面に応じて粒子に加えられる力である。
【００２１】
以上の４つの力を考慮して、上記ステップＳ６で、各粒子に加わるこれらの４つの力を全て計算した後に、次の時刻の各粒子の位置及び速度の計算を行い、計算結果に基づき、各粒子の位置及び速度などをステップＳ７で変更する。これを以下に具体的に説明する。
まず、１ｃｃ中に含まれる液体の粒子がＮ個であると仮定し、ｉを１〜Ｎまでの整数としたとき、対象とする粒子を粒子ｉで表す。図３において、ステップＳ２１で１番目の粒子について考えるため、ｉ＝１とする。
次いで、ステップＳ２２では、粒子ｉと該粒子ｉの近傍の粒子との中心間の距離ｒを計算する。すなわち、図４（Ａ）に示すように、粒子ｉと粒子ｊ、粒子ｉと粒子ｊ＋１、……、粒子ｉと粒子Ｎの各中心間距離を計算する。
【００２２】
次に、上記求められた中心間距離に基づき、ステップＳ２３において、図４（Ｂ）に示すように、粒子間における相互作用を考慮するため、各粒子から受ける力Ｆ_１ｉの計算を「数１」により行う。ここで、粒子ｉが粒子ｊから受ける力をＦｉｊで表し、逆に、粒子ｊが粒子ｉから受ける力をＦｊｉで表す。
【数１】

次に、ステップＳ２４では、各粒子ｉの中心から該粒子ｉの近傍の固体表面までの最短距離を計算する。例えば、図４（Ｃ）に示すように、粒子ｉの中心と固体表面Ｌとの最短距離ｒ_ｉＬ、粒子ｉと固体表面Ｌ＋１との最短距離ｒ_ｉＬ＋１を求める。
【００２３】
次に、上記求められた最短距離に基づき、ステップＳ２５において、各粒子ｉが固体表面から受ける力Ｆ_２ｉを計算する。例えば、図４（Ｄ）に示すように、粒子ｉが固体表面Ｌから受ける力Ｆｉ_Ｌ、粒子ｉが固体表面Ｌ＋１から受ける力Ｆｉ_Ｌ＋１などを求める。
次に、ステップＳ２６では、各粒子ｉが重力により受ける力Ｆ_３ｉ（＝ｍｇ）を計算する。ここで、１粒子の質量をｍｇとする。
さらに、ステップＳ２７では、各粒子ｉがその他の外力より受ける力Ｆ_４ｉを計算する。
次に、ステップＳ２８では、各粒子ｉが受ける上記求められた４つの力の総和（Ｆｉ＝Ｆ_１ｉ＋Ｆ_２ｉ＋Ｆ_３ｉ＋Ｆ_４ｉ）を計算する。
次に、ステップＳ２９では、ｉをｉ＋１として、次の粒子を対象として考える。すなわち、ステップＳ２９の後、ステップＳ３０でｉ＜Ｎか否か、言い換えれば、総ての粒子についての計算が終了したか否か判断し、終了していないならば、ステップＳ２２に戻り、その粒子についてステップＳ２２〜Ｓ２８の計算を行う。終了しているならば、ステップＳ３１で、各粒子のΔｔ秒後の速度ｖｉを計算する。すなわち、現在の速度をｖｉ_０とすると、ｖｉ＝ｖｉ_０＋（Ｆｉ／ｍ）×Δｔを計算する。
次に、ステップＳ３２で各粒子のΔｔ秒後の位置を計算する。以上により、ステップＳ６の計算を終了し、その計算結果に基づき、図１のステップＳ７で各粒子の位置を変更する。
【００２４】
一方、図６（Ａ），（Ｂ）は、粒子の密度計算を行う際の考え方を示している。密度は、ある領域（体積Ｖ）内に存在する存在粒子数（Ｎ）により決定される。よって、ステップＳ５０において、ある領域（Ｖ）内の粒子２０の粒子数（Ｎ）を計算する。一方、圧力は容器中に平板を置き、その平板に加わる粒子の力の総和を計算し、平板面積で割られたものを圧力とする。よって、ステップＳ５１において、上記領域（Ｖ）内のＮ個の粒子２０の質量ｍの合計（ｍ×Ｎ）を計算する。この結果、ステップＳ５２において、密度ρは、ρ＝（ｍＮ／Ｖ）により求めることができる。
また、各粒子の温度は、内部エネルギーを温度エネルギーに変換し、外から受ける温度及び外へ逃げる温度の計算は、従来通りの熱伝導計算、又は、熱伝達計算と同様な計算方法で算出することができる。
【００２５】
図７，８は、液体を貯える同一容器から同一初期条件で、粘度の異なる別々の液体が上記容器から流れ出るときの違いを示したものである。図７が粘度の低い液体を流した図であり、図８が粘度の高い液体を流した図である。図７，８において、８１は液滴、８２は液体を貯える容器、８３は床面を示す。液体が吐出される際に、容器８２の出口部に液体が付着して液体が溜まる。その後、液体の表面張力が液体の自重に耐えられなくなり、下方の床面に液体が落下していく。このとき、液体の各粒子に加わる力Ｆは、「数２」で表現される。
【数２】
Ｆ＝Ｆ_１＋Ｆ_２＋Ｆ_３
ここで、Ｆ_１は重力、Ｆ_２は粒子間力、Ｆ_３は固体表面から受ける力である。
【００２６】
粒子ｉと粒子ｊ間の相互作用による力をｆ_ｉｊとすると、「数３」のように表される。
【数３】
ｆ_ｉｊ＝ｆ（ｒ_ｉｊ）
ここで、ｆ_ｉｊは、粒子間距離ｒ_ｉｊの関数として表現される。
例えば、図１０に示すように、ある粒子間に働く曲線が２次式からなる関数で表現されるとすると、「数４」のように表される。
【数４】
ｒ_ｉｊ＜ｒ_１のときｆ_ｉｊ＝ −ｕ（ｒ_ｉｊ−ｒ_１）（ｒ_ｉｊ−ｒ_２）
ｒ_１≦ｒ_ｉｊ＜ｒ_２のときｆ_ｉｊ＝ｕ（ｒ_ｉｊ−ｒ_１）（ｒ_ｉｊ−ｒ_２）
ｒ_２≦ｒ_ｉｊのときｆ_ｉｊ＝０
ここで、ｆは力の次元、ｒは距離の次元、ｕは距離の次元を掛けることにより力の次元となるような次元を持つパラメータ（係数）である。
【００２７】
このとき、上記したように、粒子間距離が非常に小さくなると（ｆ_ｉｊ＜０のとき）、粒子間には、互いに反発する斥力が働き、一定の距離が保たれると（ｆ_ｉｊ＞０）お互いに引き合う引力が働く。粒子間の距離ｒが、横軸ｒ_ｉｊと交差する点の臨海点距離ｒ_１と等しくなると粒子間に働く力の大きさが零となり、安定した状態となる。このとき、臨海点距離ｒ_１の大きさを変更すると液体全体の密度や粘度が変更され、異なる物性を示す液体を表現することが可能となる。臨海点距離ｒ_１と距離ｒ_２の両方を変更したり、パラメータｕもこれに加えて変更することにより、種々のモデルに対応させることが可能である。
【００２８】
図９に、上記パラメータｕを変更したときの液体の状態の変化を示している。例えば、距離ｒ_２を固定した状態（例えば距離ｒ_２を３．６に固定した状態）で、粒子間の引力に関係するパラメータｕ又は臨海点距離ｒ_１又はその両方を変更すると、２次式の関数で表される曲線の高さが大きくなったり、小さくなったりして、高粘度の液体から低粘度の液体の挙動を示すことが可能となり、さらには、固体から気体の性質までを表現することが可能となる。例えば、パラメータｕを大きくすれば、粒子間又は粒子と固体表面との間での引力が大きくなり、粘度が大きくなる。一方、パラメータｕを小さくすれば、粒子間又は粒子と固体表面との間での引力が小さくなり、粘度が小さくなる。具体的には、図９において、中央のＣの粘土状態において、臨海点距離ｒ_１を大きくすれば、Ｃの粘土状態からＦの状態のように「ムース状」の粘土となる。逆に、臨海点距離ｒ_１を小さくすれば、Ｃの状態からＧの状態のように「縮こまる」状態の粘土となる。また、パラメータｕを小さくすれば、Ｃの粘土状態からＢの液体状態に変化し、さらに、パラメータｕを小さくすれば、Ａの気体状態となる。一方、逆に、パラメータｕを大きくすれば、Ｃの粘土状態からさらに粘度が大きくなったＤの粘土状態となり、途中で固まって落下しないようになる。さらにパラメータｕを大きくすれば、最初に小さく縮まってしまい、小さな球になった後にまとまって落下するようなＥの固体状態となる。
【００２９】
同様に、粒子が固体表面から、受ける力Ｆ_３においても、液と固体表面との相性や固体表面の粗さ及び加工精度の違いにおいて、粒子と固体表面との距離ｒ_ｉｋで親水性や撥水性の表面形状を表現することが可能となる。
粒子ｉと固体表面ｋ間の相互作用による力をｆ_ｉｋとすると、「数５」のように表される。
【数５】
ｆ_ｉｋ＝ｆ（ｒ_ｉｋ）
すなわち、ｆ_ｉｋは粒子と固体表面間距離ｒ_ｉｋの関数として表現される。
例えば、図１１に示すように、ある粒子と固体表面間に働く曲線を２次式からなる関数で表現されるとすると、「数６」のように表される。
【数６】
ｒ_ｉｋ＜ｒ_１のときｆ_ｉｋ＝−ｕ（ｒ_ｉｋ−ｒ_１）（ｒ_ｉｋ−ｒ_２）
ｒ_１≦ｒ_ｉｋ＜ｒ_２のときｆ_ｉｋ＝ｕ（ｒ_ｉｋ−ｒ_１）（ｒ_ｉｋ−ｒ_２）
ｒ_２≦ｒ_ｉｋのときｆ_ｉｋ＝０
このとき、ｆ_ｉｋ＞０となり粒子と固体表面との引き合う引力が強くなると粒子が表面に付着し易くなり、ｆ_ｉｋ＜０となり弱くなると斥力が発生してハジキ易くなる。
【００３０】
このように、本実施形態によれば、粒子間に働く相互作用を粒子間の距離に関するパラメータｕで、液体の物性である粘度や表面張力や密度等の値を表現することが可能となる。
以下に、具体例を説明する。
図１２は、ある液滴を固体表面上に落とした際の液体の表面形状を表現している。図１２の左側と右側では、同一物性の液体を使用しているが、粒子と固体表面間に働く力の関数を変更することで、左側に向かうに従いヌレ易い状態を示し、右側に向かいに従いヌレにくい状態を表現している。これは、同一の液体においてもガラス表面に液滴を落としたときと、ポリエチレン系の素材に落としたときの違いを表現しており、前者はヌレ易く、後者はヌレにくいことになる。
一方、図１３においては、図１２の液体が経時的変化で、ｔ＝１の状態からｔ＝２０の最終状態になるまでの変化を示しており、時間的変化が液体の粘度と関係していることを示すものである。
【００３１】
さらに、図１４は、様々な液の流速変化等に伴う液の粘度の変化を表した表である。同一の液体であっても液体の内部温度や液体の流速によって、粘度が変化する性質は、一般に、よく知られている。例えば、図１４に示すように、粘度ηとせん断速度との関係において、せん断速度にかかわらず粘度が一定なものはニュートン性があり、せん断速度が大きくなるに従い粘度が小さくなるものは非ニュートン性でかつチクソトロピーを示し、粘度が大きくなるものは非ニュートン性でかつレオペキシーを示す。これと同様なことを再現するために、あらかじめ液体の種類により、流速変化や温度変化や時間変化に伴い異なる粒子間距離の関数を持ち、各粒子の流速変化や温度変化や時間的変化を観測制御することにより、通常と異なる粒子間の相互作用を行う。これにより、時間的変化に伴い粘度が変更し、実行したものが、図１５，１６である。
【００３２】
図１５において１４０、１４１、１４２、図１６の１４３、１４４、１４５の順に進むほど、時間経過に伴い液体の粘度が変化しており、図１５の１４０及び１４１の状態では、液体が高粘度状態であり、液体が液体容器の固体表面に付着していることを示す。１４２及び１４３の状態では、液体が比較的低粘度状態であり、液体が固体表面に付着しておらず、順次落下することを示す。さらに、図１６の１４４及び１４５の状態では、液体が再度高粘度状態となり、液体が固体表面に付着していることを示す。
さらに、図１７は、毛細管現象を表しており、液体を溜める容器内に細い管を挿入したとき、容器内の液体がその細い管を伝わって毛細管現象により上昇する状態を示している。
また、図１８は、あらかじめ液体を固体表面上に落として置き、その後、落とす位置をずらせて、さらに液体を固体表面上に落とした状態を示している。液体が、合体及び分割の状態を頻繁に起こす液体においては、従来の流体解析手法では計算不可能であった、界面の状態変化を計算することが容易となる。
【００３３】
図２３は、粒子数４００個が１ｃｃとしたときのある液体の高さ比（約１ｃｃの液体を数滴落下させたときの各液滴の高さのうちの最大高さで各液滴の高さを割った値）又は直径の比（各液滴の直径のうち最大直径で各液滴の直径を割った値）と体積比（各液滴の体積のうちの最大体積で各液滴の体積を割った値）との関係を示すグラフである。図２４はその結果を示す説明図である。図２３中、直径比と体積比について、一点鎖線Ｉは上記実施形態によるシミュレーションしたときの関係を示すグラフ、実線ＩＩは実験により測定して求められた関係を示すグラフである。また、高さ比と体積比にいて、太線ＩＩＩは上記実施形態によるシミュレーションしたときの関係を示すグラフ、実線ＩＩは実験により測定して求められた関係を示すグラフである。図２４において、（Ａ）は粒子数が４００個で液体約１ｃｃに相当する場合であり、実験により測定する実際の約１ｃｃの液体の状態を（Ｂ）に示す。同様に、（Ｃ）は粒子数が２００個で液体約０．５ｃｃに相当する場合であり、実験により測定する実際の約０．５ｃｃの液体の状態を（Ｄ）に示す。同様に、（Ｅ）は粒子数が４０個で液体約０．１ｃｃに相当する場合であり、実験により測定する実際の約０．１ｃｃの液体の状態を（Ｆ）に示す。図２３，２４より、実験により求められた関係とシミュレーションによる関係とが近似していることがわかる。
【図面の簡単な説明】
【図１】本発明の一実施形態にかかる粒子型流体シミュレーション方法を示すフローチャートである。
【図２】上記シミュレーション方法を実施するためのシミュレーション装置のブロック図である。
【図３】上記実施形態におけるより具体的な粒子の位置及び速度の計算アルゴリズムを示す図である。
【図４】（Ａ）〜（Ｄ）はそれぞれ図３のアルゴリズム中の距離計算動作の説明図である。
【図５】上記実施形態において、流体の不安定な分子間力を粒子間同士の相互作用で近似する考え方を示す図である。
【図６】（Ａ）及び（Ｂ）は上記実施形態における密度の計算方法を示すフローチャート及びその説明図である。
【図７】上記実施形態における液体のある物性のときの液体の落下状況の比較図である。
【図８】上記実施形態における液体の物性が図７と異なるときの液体の落下状況の比較図である。
【図９】上記実施形態における粒子間の距離に依存したパラメータ値又は距離ｒ_１を変更させたときの液体の落下状況の比較図である。
【図１０】上記実施形態における粒子間の相互作用を決定する粒子間距離に依存した関数図である。
【図１１】上記実施形態における粒子と固体表面との相互作用を決定する粒子と固体表面の最短距離に依存した関数図である。
【図１２】上記実施形態における粒子と固体表面との距離に依存したパラメータを変更したときの比較図である。
【図１３】上記実施形態における液体形状の経時的変化図である。
【図１４】上記実施形態における様々な液体の速度変化、時間的変化に伴う液粘度の変化グラフである。
【図１５】上記実施形態における時間的変化に伴い粘度が変更した場合の液の落下図である。
【図１６】図１５に続く、上記実施形態における時間的変化に伴い粘度が変更した場合の液体の落下図である。
【図１７】上記実施形態における毛細管現象図である。
【図１８】上記実施形態における液体の合体図である。
【図１９】（Ａ），（Ｂ）は容器中に入れた液体の容器を傾ける前及び傾けたときの液体界面の経時的変化図である。
【図２０】管内の液体の流れ図である。
【図２１】（Ａ）〜（Ｄ）はＶＯＦ法による計算方法を示し、（Ａ）は液体が容器内に入った状態の斜視図、（Ｂ）は（Ａ）の一部拡大図、（Ｃ），（Ｄ）は同一占有率の比較図である。
【図２２】（Ａ），（Ｂ）は流体の流れと液のヌレによる横飛び現象が生じているときの説明図及びその拡大図である。
【図２３】ある液体を数滴落下させたときの液体の体積比と直径比又は高さ比との関係について、上記実施形態にかかるシミュレーションにより求めた関係と実験により測定した関係とを比較するグラフである。
【図２４】図２３においてシミュレーションにおける粒子数４００、２００、４０個の場合の状態とそれぞれの実際の液体の状態を示す比較図である。
【符号の説明】
１入力部
２記憶部
３演算部
４出力部
８１液滴
８２液を貯える容器
８３床面
１４０、１４１高粘度状態で液体が固体表面に付着している状態
１４２、１４３比較的低粘度の液体に変更したときの付着状態
１４４、１４５再度高粘度の液体に変更したときの付着状態[0001]
TECHNICAL FIELD OF THE INVENTION
The present invention provides a particle-type fluid simulation method for quantitatively expressing a change in the shape of a fluid (liquid or gas) or a surface shape (interface state) from a fluid to a solid, and predicting dynamic and static fluid behavior. And an apparatus.
[0002]
[Prior art]
Conventional fluid analysis techniques can be classified into two types of methods. One is to solve Navier-Stokes equation (Navier-Stokes equation) which is a basic equation of fluid, and the other is to solve a free surface. This is a VOF method (A Volume of Fluid Algorithm) that is handled. The former does not mention the surface shape of the fluid, it calculates the flow velocity, pressure and fluid resistance inside the fluid, and the latter can calculate the surface shape of the fluid, but does not calculate it accurately There is such a problem.
[0003]
A specific description will be given with reference to FIGS. First, in the former method, even if the liquid itself changes with time, such change is not calculated. As shown in FIGS. 19A and 19B, when the container 60 containing the liquid 62 is inclined obliquely, the shape of the whole liquid changes according to the inclination of the container 60, but this kind of problem occurs. There is a problem that cannot be dealt with. In the figure, 63 is the interface of the liquid 62 in the container 60. However, as shown in FIG. 20, calculation can be performed to solve the velocity v, the flow trace, the pressure difference ΔP, and the fluid resistance of the liquid flowing through the entire pipe. The surface shape itself such as coalescence cannot be handled. The latter method will be described with reference to FIG. 21. In the method of calculating the shape of the liquid 72 based on the volume occupancy in the divided area as shown in FIG. Regarding the boundary of the liquid 72 in 70 (the interface 71 of the liquid 72), as long as the volume occupancy is the same in the divided areas as shown in FIGS. There is a problem that cannot be determined exactly. Therefore, it is not possible to accurately calculate the surface shape of the liquid 72. Similarly, surface shape predictions such as liquid separation and coalescence cannot be addressed.
[0004]
However, in recent years, a technique for forming a film accurately and uniformly on a solid surface using various liquids has become indispensable not only for semiconductors but also for various technical fields. At this time, predicting what the surface shape of the liquid will be has become very important. For example, as shown in FIG. 22A, not only the entire flow velocity v of the liquid 92 discharged from the container 91 by the action of the pressure P, but also the relationship between the liquid surface and the solid surface of the container 91 due to liquid wetting or repelling. The phenomenon in which the liquid 92 causes the horizontal jump 94 has a problem that it is difficult to predict with the conventional analysis method for the above-described reason. In FIG. 22B, reference numeral 93 denotes a state in which liquid droplets adhere to the opening end surface of the container 91 due to liquid wetting.
[0005]
[Problems to be solved by the invention]
In the conventional fluid analysis technology and free surface analysis technology typified by the VOF method, as described in the above-described specific example, the calculation of the micro droplet shape itself and the swelling and cissing of the liquid and solid surface shapes are involved. The phenomenon that the liquid jumped sideways was not able to be calculated. Therefore, the only way to reproduce such a phenomenon was to repeat trial and error experiments.
Accordingly, an object of the present invention is to solve the above-mentioned problems of the conventional fluid analysis technology and the free surface analysis technology represented by the VOF method, and to calculate the microscopic surface shape itself of the fluid. It is an object of the present invention to provide a particle-type fluid simulation method and apparatus capable of more accurately calculating a phenomenon associated with swelling, swelling and cissing.
[0006]
Means for Solving the Problems and Their Effects
In order to achieve the above object, the present invention is configured as follows.
That is, according to the particle-type fluid simulation method according to the first aspect of the present invention, the fluid is not a continuum, but each minute region is a virtual particle.AndThe interaction of forces acting on each particleExpressBased on a function that uses the distance between particles as a parameter, the change over time in the position of each particle is determined, and the change in the fluid shape is reproduced.Thing,
The above function comprises a curve or a straight line characterized by generating repulsive forces when the distance between the respective particles becomes a certain critical point distance and generating an attractive force when the distance becomes a certain distance beyond the critical point distance, This curve or straight line changes the temperature, velocity, or time of each particle, and the curve or straight line itself changes, thereby expressing a change in non-Newtonian viscosity and reproducing a change in the shape of the fluid.Like that.
In the invention according to the second aspect of the present invention, in the first aspect,SaidBased on the function,eachIt can also determine changes in particle pressure over time and reproduce changes in fluid shape.it can.
[0007]
The present invention3According to the particle type fluid simulation method according to the aspect, the fluid is not a continuum, but each minute region is a virtual particle.ageIn order to reproduce hydrophilicity and water repellency in the relationship between each particle and the solid surface, the interaction between the force of each particle and the solid surface is determined based on a function using the distance between the solid surface and the particle as a parameter. The change with time of the position of the particles is determined to reproduce phenomena such as hydrophilicity and water repellency between the fluid and the solid surface.
The present invention4According to the aspect, the second3In embodiments,SaidBased on the function,eachThe change with time of the pressure of the particles can also be determined, and phenomena such as hydrophilicity and water repellency between the fluid and the solid surface can be reproduced.
The present invention5According to the aspect, the second3 or 4In embodiments,SaidThe function is a curve characterized by generating repulsive forces when the distance between each particle and the solid surface reaches a certain critical point distance, and generating attractive forces when a certain distance exceeds the critical point distance. It consists of a straight line, and this curve or straight line changes with the temperature change, speed change, or time change of each particle, by changing the curve or straight line itself,Non-Newtonian viscosityExpress the change ofAnd reproduces phenomena such as hydrophilicity and water repellency between the fluid and the solid surfaceIt can also be done.
[0010]
Therefore, according to the various aspects of the present invention, the fluid is expressed by replacing the fluid with virtual particles, and the interaction of forces acting between the particles is determined based on a function using the distance between the particles as a parameter. Changes in the position, velocity, density, pressure, or temperature of the particles over time can be determined.
According to the various aspects of the present invention, the interaction between the force of the particles and the solid surface is based on a function using the distance between the solid surface and the particles as a parameter, and the position, velocity, and density of each particle. By determining the change over time in pressure, temperature, and temperature, wetting (hydrophilicity) and cissing (water repellency) between the fluid and the solid surface can be reproduced.
Also, by controlling the interaction of the force acting between each particle with a temperature change, a speed change, or a function with the distance between each particle as a parameter along with the time change, the thixotropy of the non-Newtonian viscosity, etc. It is possible to cope with physical properties in which the viscosity changes with time.
Further, by having a plurality of parameters such as mass, volume, temperature, viscosity, and surface tension of each particle, a fluid having two or more different physical properties or a two-phase flow (solid and liquid, solid and gas, liquid And gas, solid, liquid, and gas) or more.
Therefore, according to the present invention, the surface shape from a fluid (gas or liquid) to a solid is calculated by the inter-particle force using a particle type model, so that the wettability due to the viscosity and surface tension of the fluid, the microscopic surface shape, and the like are obtained. It can simulate more accurately, quantitatively express the change in the shape of a fluid (liquid or gas), or the surface shape (interface state) from fluid to solid, and predict dynamic and static fluid behavior be able to. Therefore, for example, it is possible to determine the optimal shape of the fluid storage container or the ejection nozzle or the like at the time of fluid ejection such as a liquid application device and the optimal viscosity of the ejected fluid.
[0011]
BEST MODE FOR CARRYING OUT THE INVENTION
Hereinafter, embodiments according to the present invention will be described in detail with reference to the drawings.
The present application discloses a patent application No. Hei 8-178894 (named "liquid and wetting evaluation device"), which has already been filed by the present applicant, is not intended to transfer a fluid (gas or liquid) to a continuum but to each minute area. It describes a method of treating particles as virtual particles, and does not describe in detail a specific calculation method. Thus, as a result of earnest research, the present inventor has created an invention that can clarify a specific calculation algorithm when a fluid is represented by particles. The content is roughly composed of two points. The first is to determine the interaction of the forces acting between the particles based on a function using the distance between the particles as a parameter, and the second is to determine the interaction between the particles and the solid surface. This paper describes a method for reproducing the hydrophilicity and water repellency between the liquid and the solid surface by determining the interaction based on a function using the shortest distance between the particle and the solid surface as a parameter.
[0012]
In other words, in the present embodiment, the fluid is not a continuum, but each micro-region is treated as a virtual particle, and the interaction between forces acting on each of these particles is defined as a distance between the particles as a parameter. Based on the function, the change over time in the position of each particle (other than that, for example, velocity, density, pressure, temperature, etc.) is determined, and the change in the shape of the fluid is reproduced. Then, based on the function of the interaction between the forces acting on the particles and the distance between the particles as a parameter, the change with time of the pressure of the particles can also be determined, and the change in the shape of the fluid can be reproduced. Further, the function of the parameter of the distance expressing the interaction of the forces acting between the respective particles, reciprocal force occurs when the distance between the particles reaches a certain critical point distance, a constant beyond this critical point distance It consists of a curve or a straight line characterized by causing an attractive force with each other when the distance is reached, and this curve or the straight line changes with the temperature change, the speed change, or the time change of each particle. It is also possible to express the change in viscosity with time due to the thixotropy of non-Newtonian viscosity.
[0013]
Further, in the present embodiment, the fluid is not a continuum, but each microregion is treated as a virtual particle, and in order to reproduce hydrophilicity and water repellency in the relationship between each particle and the solid surface, each particle is Time-dependent changes in the position of each particle (for example, velocity, density, pressure, temperature, etc.) based on the interaction of the force between the particle and the solid surface based on a function using the distance between the solid surface and the particle as a parameter Is determined, and phenomena such as hydrophilicity and water repellency between the fluid and the solid surface are reproduced. Then, based on a function using the distance between the solid surface and the particle as a parameter, the interaction of the force between each particle and the solid surface is also determined as to the time-dependent change in the pressure of each of the particles. Phenomena such as hydrophilicity and water repellency can also be reproduced. Further, the distance parameter function expressing the interaction of the force acting between each particle and the solid surface generates repulsive forces when the distance between each particle and the solid surface reaches a certain critical point distance. A curve or a straight line characterized by generating an attractive force with each other at a certain distance beyond the critical point distance, and the curve or the straight line is accompanied by a temperature change, a speed change, or a time change of each particle. By changing the curve or the straight line itself, it is also possible to express the change in viscosity with time due to the thixotropic property of the non-Newtonian viscosity.
Each of the particles may have a plurality of particle parameters such as mass, volume, viscosity, and surface tension of each particle, so that various models can be expressed.
[0014]
Hereinafter, a particle-type fluid simulation method and apparatus according to an embodiment of the present invention will be described with reference to the drawings.
FIG. 1 is a flowchart for carrying out the particle-type fluid simulation method according to the embodiment of the present invention. FIG. 2 is a block diagram of the simulation apparatus for performing the particle fluid simulation method according to the embodiment of the present invention.
In FIG. 2, a fluid to be analyzed and analysis conditions such as boundary conditions are input from an input unit 1 such as a keyboard or a mouse and stored in the storage unit 2, and then stored in the storage unit 2 by the calculation unit 3. The simulation method is performed based on the data, and the simulation result is output by the output unit 4 such as a display screen or a printing device.
[0015]
The simulation method performed by the arithmetic unit 3 will be described with reference to the flowchart of FIG.
First, in step S1, the physical property values of the fluid or the physical property values of the fluid and the solid surface are set and input from the input unit 1 to the storage unit 2 and stored, and the position and shape of the fluid are set as the initial arrangement. Similarly, the position and shape of the fluid and solid surfaces are input from the input unit 1 and stored in the storage unit 2.
Next, in step S2, the arithmetic unit 3 searches for a neighboring particle of the particle of interest, and if there is a neighboring particle, calculates the distance between the center of the particle of interest and the neighboring particle. A nearby solid surface of the particle to be searched is searched, and if a nearby solid surface exists, the distance between the particle of interest and the nearby solid surface is calculated.
[0016]
Next, in step S3, the input unit 1 selects a parameter depending on the distance between the particles or the distance between the solid surfaces of the particles based on the calculation result of the calculation unit 3. When selecting a parameter depending on the distance between the particles, in step S4, reference is made to an interaction parameter depending on the distance between the particles according to the type of liquid or the like. For example, depending on each of the objects to be analyzed, such as between liquid A and liquid A, between liquid B and liquid B, or between liquid A and liquid B, the parameters associated with the temperature change, the speed change, Appropriate parameters are selected from the accompanying parameters and the parameters accompanying the time change. Thereafter, the process proceeds to step S6. On the other hand, when a parameter depending on the distance between the particle and the solid surface is selected in step S3 based on the calculation result in the calculation unit 3, the process proceeds to step S5. In this step S5, an interaction parameter depending on the distance between the particle and the solid surface is referred to according to the solid surface. For example, an appropriate parameter is selected from a parameter associated with a temperature change, a parameter associated with a speed change, a parameter associated with a time change, and the like in accordance with each object to be analyzed such as the solid surface C or the solid surface D. Thereafter, the process proceeds to step S6. In addition, various parameters of the interaction depending on the distance between the particles according to the type of the liquid, and various parameters between the particles and the solid surface are stored in the storage unit 2 in advance, and are appropriately stored. The user selects and refers to the data of these parameters in the storage unit 2 by an input operation from the input unit 1.
[0017]
Next, in step S6, the calculation unit 3 calculates the force received by each particle using the selected parameters.
Next, based on the calculation result in step S6, in step S7, the calculation unit 3 changes the position of each particle. At this time, the pressure, speed, density, and the like are changed as necessary, other than the position of each particle or without changing the position.
Next, in step S8, the arithmetic unit 3 determines whether or not the analysis is completed. If not, the process returns to step S2. If completed, the result is appropriately output by the output unit 4 in step S9, and the analysis is completed.
[0018]
Here, the calculation of the force applied to each particle in step S6 and the change in the position of each particle in step S7 will be described in more detail with reference to FIGS.
FIGS. 3 and 4 describe the algorithms performed in the operations of steps S6 and S7 according to the embodiment. Here, it is considered that each particle receives roughly the following four forces, and changes the position, velocity, pressure, temperature, and density of each particle at the next time compared to the current time. .
[0019]
First, the first force is a force acting between the particles, and this force makes it possible to express the viscosity and ductility inside the fluid. In the calculation method, when considering particles i and particles j, it is considered that the force F generated between the particles changes depending on the distance r between the centers of the particles i and j. In general, as shown in FIG. 5, in view of the fact that an intermolecular force acts between two molecules, it can be approximately approximated that an interaction using the distance between the particles acts as a parameter between the particles. Is what you think. Specifically, the vertical axis represents the force F generated between the particles, and the horizontal axis represents the distance r between the particles._xFIG. 5 is a graph showing the relationship between the two cases. The point at which this graph crosses the horizontal axis is the critical point distance r described below._cIt is. Therefore, in this graph, as shown in FIG. 5A, the distance r between particles, which is the distance between the centers of particles i and j, is a certain critical point distance r._cWhen the distance is shorter, repulsive forces are generated, and as shown in FIG._cWhen the force becomes equal to, no force acts between the particles, and as shown in FIG._cIf the length is longer, attractive force is generated between the particles. Further, although not specifically shown in FIG. 5, it can be seen that when the distance between the particles is further increased, no force is exerted between the particles. Also, the viscosity basically increases as the area representing the attractive force increases, and when the viscosity becomes extremely strong, the solidified physical properties are obtained. Therefore, the viscosity can be changed by changing the area representing the attractive force. Also, the critical point distance r_cBy changing, the density of the particles also changes. That is, the critical point distance r_cAs the distance becomes larger, the distance between particles becomes longer, so that the density of the particles becomes lower. This makes it possible to change the density of the particles by changing the distance between the particles.
[0020]
The second force is a force acting between the particles and the solid surface, and it is possible to express the hydrophilicity and water repellency of the liquid by this force. The calculation method calculates the shortest distance between the center of the particle and the solid surface in the same manner as the force between the particles, and expresses the distance as a function using the distance as a parameter. The function is basically the same as the calculation between particles, and the function changes depending on the compatibility between the particles and the solid surface.
The third force is gravitational force, and gravitational force always acts on all particles in a fixed direction.
The fourth force is a force applied to the particles depending on the situation, such as an electrostatic force or an inertial force, as other external forces.
[0021]
In consideration of the above four forces, in step S6, after calculating all these four forces applied to each particle, the position and velocity of each particle at the next time are calculated, and based on the calculation result, The position and velocity of each particle are changed in step S7. This will be specifically described below.
First, it is assumed that the number of liquid particles contained in 1 cc is N, and when i is an integer from 1 to N, the target particle is represented by particle i. In FIG. 3, since the first particle is considered in step S21, i = 1.
Next, in step S22, the distance r between the centers of the particle i and the particle near the particle i is calculated. That is, as shown in FIG. 4A, the center distances of the particle i and the particle j, the particle i and the particle j + 1,..., The particle i and the particle N are calculated.
[0022]
Next, based on the calculated center-to-center distance, in step S23, as shown in FIG. 4B, in order to consider the interaction between the particles, the force F received from each particle is considered.₁The calculation of i is performed by “Equation 1”. Here, the force that particle i receives from particle j is represented by Fij, and the force that particle j receives from particle i is represented by Fji.
(Equation 1)

Next, in step S24, the shortest distance from the center of each particle i to the solid surface near the particle i is calculated. For example, as shown in FIG. 4C, the shortest distance r between the center of the particle i and the solid surface L_iL, The shortest distance r between the particle i and the solid surface L + 1_{iL + 1}Ask for.
[0023]
Next, based on the shortest distance obtained above, in step S25, the force F received by each particle i from the solid surface is obtained.₂Calculate i. For example, as shown in FIG. 4D, a force Fi that a particle i receives from a solid surface L_L, The force Fi that the particle i receives from the solid surface L + 1_{L + 1}Ask for.
Next, in step S26, the force F applied to each particle i by gravity is obtained.₃Calculate i (= mg). Here, the mass of one particle is mg.
Further, in step S27, the force F that each particle i receives from other external force₄Calculate i.
Next, in step S28, the sum of the four forces obtained above that each particle i receives (Fi = F₁i + F₂i + F₃i + F₄Calculate i).
Next, in step S29, i is set to i + 1, and the next particle is considered. That is, after step S29, it is determined whether or not i <N in step S30, in other words, whether or not the calculation has been completed for all the particles. If not, the process returns to step S22 and returns to step S22. Are calculated in steps S22 to S28. If it has ended, in step S31, the velocity vi of each particle after Δt seconds is calculated. That is, the current speed is set to vi₀Then, vi = vi₀Calculate + (Fi / m) × Δt.
Next, in step S32, the position of each particle after Δt seconds is calculated. Thus, the calculation in step S6 is completed, and the position of each particle is changed in step S7 in FIG. 1 based on the calculation result.
[0024]
On the other hand, FIGS. 6A and 6B show the concept when calculating the density of particles. The density is determined by the number (N) of existing particles existing in a certain region (volume V). Therefore, in step S50, the number (N) of the particles 20 in a certain region (V) is calculated. On the other hand, for the pressure, a flat plate is placed in a container, the total force of the particles applied to the flat plate is calculated, and the pressure divided by the flat plate area is defined as the pressure. Therefore, in step S51, the total (m × N) of the masses m of the N particles 20 in the region (V) is calculated. As a result, in step S52, the density ρ can be obtained from ρ = (mN / V).
In addition, the temperature of each particle converts internal energy into temperature energy, and the calculation of the temperature received from the outside and the temperature escaping outside is calculated by a conventional heat conduction calculation, or a calculation method similar to the heat transfer calculation. be able to.
[0025]
7 and 8 show the difference when different liquids having different viscosities flow out of the same container under the same initial conditions under the same initial conditions. 7 is a diagram in which a low-viscosity liquid flows, and FIG. 8 is a diagram in which a high-viscosity liquid flows. 7 and 8, reference numeral 81 denotes a droplet, 82 denotes a container for storing a liquid, and 83 denotes a floor surface. When the liquid is discharged, the liquid adheres to the outlet of the container 82 and accumulates. Thereafter, the surface tension of the liquid cannot withstand the weight of the liquid, and the liquid falls on the floor below. At this time, the force F applied to each particle of the liquid is expressed by “Equation 2”.
(Equation 2)
F = F₁+ F₂+ F₃
Where F₁Is gravity, F₂Is the interparticle force, F₃Is the force received from the solid surface.
[0026]
Let the force due to the interaction between particle i and particle j be f_ijThen, it is expressed as “Equation 3”.
(Equation 3)
f_ij  = F (r_ij)
Where f_ijIs the distance r between particles_ijExpressed as a function of
For example, as shown in FIG. 10, if a curve acting between certain particles is represented by a function consisting of a quadratic equation, it is expressed as “Equation 4”.
(Equation 4)
r_ij<R₁  When f_ij  = -U (r_ij-R₁) (R_ij-R₂)
r₁≤r_ij<R₂  When f_ij  = U (r_ij-R₁) (R_ij-R₂)
r₂≤r_ij  When f_ij  = 0
Here, f is the dimension of force, r is the dimension of distance, and u is a parameter (coefficient) having a dimension that becomes the dimension of force by multiplying by the dimension of distance.
[0027]
At this time, as described above, if the distance between particles becomes very small (f_ij<0), repulsive forces repelling each other act between the particles, and when a certain distance is maintained (f_ij> 0) Attraction to each other works. The distance r between particles is represented by the horizontal axis r_ijSeaside distance r at the intersection with₁When it becomes equal, the magnitude of the force acting between the particles becomes zero, and a stable state is obtained. At this time, the critical point distance r₁When the size of the liquid is changed, the density and viscosity of the entire liquid are changed, and it is possible to express liquids having different physical properties. Waterfront point distance r₁And distance r₂By changing both of them, and by changing the parameter u in addition to this, it is possible to correspond to various models.
[0028]
FIG. 9 shows a change in the state of the liquid when the parameter u is changed. For example, distance r₂Is fixed (for example, the distance r₂Is fixed to 3.6), the parameter u related to the attractive force between the particles or the critical point distance r₁If both or both are changed, the height of the curve represented by the quadratic function increases or decreases, and it becomes possible to show the behavior of a high-viscosity liquid to a low-viscosity liquid, Can express everything from solid to gaseous properties. For example, if the parameter u is increased, the attractive force between the particles or between the particles and the solid surface increases, and the viscosity increases. On the other hand, if the parameter u is reduced, the attractive force between the particles or between the particles and the solid surface decreases, and the viscosity decreases. Specifically, in FIG. 9, in the clay state of the center C, the critical point distance r₁Is increased, the clay state of C becomes a mousse-like clay like the state of F. Conversely, the critical point distance r₁Is smaller, the clay is in a state of “shrinking” from the state of C to the state of G. If the parameter u is reduced, the state changes from the clay state of C to the liquid state of B. If the parameter u is further reduced, the state changes to the gas state of A. On the other hand, if the parameter u is increased, the clay state of C becomes a clay state of D having a further increased viscosity from the clay state of C. If the parameter u is further increased, the solid state is reduced to a small sphere at first, and becomes a solid state of E such that the sphere becomes a small ball and then falls together.
[0029]
Similarly, the force F exerted by the particle from the solid surface₃In the case of the difference between the compatibility between the liquid and the solid surface, the roughness of the solid surface, and the processing accuracy, the distance r between the particles and the solid surface_ikMakes it possible to express a hydrophilic or water-repellent surface shape.
Let the force due to the interaction between particle i and solid surface k be f_ikThen, it is expressed as “Equation 5”.
(Equation 5)
f_ik  = F (r_ik)
That is, f_ikIs the distance r between the particle and the solid surface_ikExpressed as a function of
For example, as shown in FIG. 11, if a curve acting between a certain particle and a solid surface is represented by a function consisting of a quadratic equation, it is expressed as “Equation 6”.
(Equation 6)
r_ik<R₁  When f_ik= -U (r_ik-R₁) (R_ik-R₂)
r₁≤r_ik<R₂  When f_ik= U (r_ik-R₁) (R_ik-R₂)
r₂≤r_ik  When f_ik= 0
At this time, f_ik> 0 and the attractive force between the particles and the solid surface increases, the particles tend to adhere to the surface, and f_ikIf it becomes <0 and becomes weak, a repulsive force is generated and cissing becomes easy.
[0030]
As described above, according to the present embodiment, the interaction between the particles can be expressed by the parameter u relating to the distance between the particles, such as viscosity, surface tension, and density, which are physical properties of the liquid.
Hereinafter, a specific example will be described.
FIG. 12 illustrates a surface shape of a liquid when a certain droplet is dropped on a solid surface. Although liquids having the same physical properties are used on the left and right sides in FIG. 12, by changing the function of the force acting between the particles and the solid surface, the liquid tends to swell toward the left side, and slims toward the right side. Expresses a difficult state. This expresses the difference between the case where the same liquid is dropped on a glass surface and the case where the same liquid is dropped on a polyethylene-based material.
On the other hand, FIG. 13 shows the change of the liquid of FIG. 12 with time, from the state at t = 1 to the final state at t = 20, and the temporal change is related to the viscosity of the liquid. It indicates that there is.
[0031]
FIG. 14 is a table showing the change in the viscosity of the liquid accompanying the change in the flow rate of various liquids. It is generally well known that the viscosity of the same liquid changes depending on the internal temperature of the liquid and the flow velocity of the liquid. For example, as shown in FIG. 14, in the relationship between the viscosity η and the shear rate, those having a constant viscosity irrespective of the shear rate have Newtonian properties, and those having the viscosity decrease as the shear rate increases increase the non-Newtonian properties. Those exhibiting high viscosities and thixotropy and having a high viscosity are non-Newtonian and exhibit rheopexy. In order to reproduce the same thing, we have a function of the distance between particles depending on the flow velocity change, temperature change and time change depending on the type of liquid in advance, and observe the flow velocity change, temperature change and time change of each particle By controlling, unusual interaction between particles is performed. As a result, the viscosity is changed in accordance with the temporal change, and the results are shown in FIGS.
[0032]
In FIG. 15, the viscosity of the liquid changes with time as it proceeds in the order of 140, 141, 142 and 143, 144, 145 in FIG. 16. In the state of 140 and 141 in FIG. Which indicates that the liquid is attached to the solid surface of the liquid container. The states of 142 and 143 indicate that the liquid is in a relatively low viscosity state, that the liquid is not attached to the solid surface, and that the liquid falls sequentially. Further, the states of 144 and 145 in FIG. 16 indicate that the liquid is in the high viscosity state again, indicating that the liquid has adhered to the solid surface.
Further, FIG. 17 illustrates a capillary phenomenon, in which when a thin tube is inserted into a container for storing liquid, the liquid in the container rises due to the capillary phenomenon along the thin tube.
FIG. 18 shows a state in which the liquid is dropped on the solid surface in advance, and then the drop position is shifted, and the liquid is further dropped on the solid surface. In a liquid in which the liquid frequently undergoes a coalescing and splitting state, it becomes easy to calculate a state change of the interface, which cannot be calculated by a conventional fluid analysis method.
[0033]
FIG. 23 shows the height ratio of a certain liquid when 400 particles are 1 cc (the maximum height of each droplet when several drops of approximately 1 cc of liquid are dropped). The height ratio) or the ratio of the diameters (the diameter of each droplet divided by the maximum diameter of the droplets) and the volume ratio (the maximum volume of the droplets) 2 is a graph showing a relationship between a value obtained by dividing a volume of the sample. FIG. 24 is an explanatory diagram showing the result. In FIG. 23, with respect to the diameter ratio and the volume ratio, an alternate long and short dash line I is a graph showing a relationship when the simulation is performed according to the above-described embodiment, and a solid line II is a graph showing a relationship measured and obtained by an experiment. In the height ratio and the volume ratio, a thick line III is a graph showing a relationship when the simulation is performed according to the above-described embodiment, and a solid line II is a graph showing a relationship obtained by measurement through experiment. 24A shows a case where the number of particles is 400 and corresponds to about 1 cc of liquid. FIG. 24B shows an actual state of about 1 cc of liquid measured by an experiment. Similarly, (C) is a case where the number of particles is 200 and corresponds to about 0.5 cc of liquid, and the actual state of about 0.5 cc of liquid measured by experiment is shown in (D). Similarly, (E) is a case where the number of particles is 40 and corresponds to about 0.1 cc of the liquid. The actual state of the liquid of about 0.1 cc measured by an experiment is shown in (F). 23 and 24 that the relationship obtained by the experiment and the relationship obtained by the simulation are similar.
[Brief description of the drawings]
FIG. 1 is a flowchart illustrating a particle-type fluid simulation method according to an embodiment of the present invention.
FIG. 2 is a block diagram of a simulation device for performing the simulation method.
FIG. 3 is a diagram showing a more specific algorithm for calculating the position and velocity of particles in the embodiment.
FIGS. 4A to 4D are explanatory diagrams of a distance calculation operation in the algorithm of FIG. 3;
FIG. 5 is a view showing a concept of approximating an unstable intermolecular force of a fluid by an interaction between particles in the embodiment.
FIGS. 6A and 6B are a flowchart and a diagram illustrating a density calculation method in the embodiment.
FIG. 7 is a comparison diagram of a state of dropping of the liquid when the liquid has certain physical properties in the embodiment.
FIG. 8 is a comparison diagram of a dropping state of the liquid when the physical properties of the liquid are different from those of FIG. 7 in the embodiment.
FIG. 9 shows a parameter value or a distance r depending on a distance between particles in the embodiment.₁FIG. 9 is a comparison diagram of a liquid drop state when the number is changed.
FIG. 10 is a function diagram depending on an interparticle distance for determining an interaction between particles in the embodiment.
FIG. 11 is a function diagram that depends on the shortest distance between the particle and the solid surface to determine the interaction between the particle and the solid surface in the embodiment.
FIG. 12 is a comparison diagram when a parameter depending on a distance between a particle and a solid surface in the embodiment is changed.
FIG. 13 is a diagram showing a change over time in a liquid shape in the embodiment.
FIG. 14 is a graph showing a change in liquid viscosity with a change in speed and a change with time of various liquids in the embodiment.
FIG. 15 is a drop diagram of the liquid when the viscosity changes with time in the embodiment.
FIG. 16 is a continuation diagram of FIG. 15 of the liquid drop when the viscosity changes with time in the embodiment.
FIG. 17 is a diagram of a capillary phenomenon in the embodiment.
FIG. 18 is a combined view of liquids in the embodiment.
FIGS. 19A and 19B are diagrams showing the change over time of the liquid interface before and when the liquid container placed in the container is tilted.
FIG. 20 is a flow diagram of a liquid in a tube.
21A to 21D show a calculation method by the VOF method, FIG. 21A is a perspective view of a state in which a liquid is contained in a container, FIG. 21B is a partially enlarged view of FIG. (C) and (D) are comparison diagrams of the same occupancy.
FIGS. 22A and 22B are an explanatory view and an enlarged view thereof when a horizontal jump phenomenon occurs due to the flow of a fluid and the wetting of the liquid.
FIG. 23 compares the relationship between the volume ratio and the diameter ratio or the height ratio of the liquid when a few drops of the liquid are dropped, with the relationship obtained by the simulation according to the embodiment and the relationship measured by experiment. It is a graph.
FIG. 24 is a comparison diagram showing the state in the case of 400, 200, and 40 particles in the simulation in FIG. 23 and the actual state of each liquid.
[Explanation of symbols]
1 Input section
2 Storage unit
3 Operation part
4 Output section
81 droplets
82 Container for storing liquid
83 floor
140, 141 High viscosity state with liquid attached to solid surface
142, 143 Adhesion state when changed to liquid with relatively low viscosity
144, 145 Adhesion state when changed to high viscosity liquid again

Claims

流体を、連続体ではなく、各微小領域を仮想的な粒子とし、各粒子に働く力の相互作用を表現する粒子間の距離をパラメータとした関数に基づいて、各粒子の位置の経時的変化を決定し、流体の形状の変化を再現するものであって、
上記関数は、各粒子間の距離がある臨界点距離となるとお互いに斥力を生じ、この臨界点距離を越えて一定の距離になるとお互いに引力を生じることを特徴とする曲線又は直線からなり、この曲線又は直線は、各粒子の温度変化、速度変化、又は時間変化に伴い、曲線又は直線自体が変化することにより、非ニュートン性粘度の変化を表現し、流体の形状の変化を再現するようにした粒子型流体シミュレーション方法。A fluid, rather than a continuum, with reference to each micro area as virtual particles, a function that the distance was the parameter between the particles representing the interaction force acting on each particle, with time the position of each particle Determine the change and reproduce the change in the shape of the fluid ,
The above function comprises a curve or a straight line characterized by generating repulsive forces when the distance between the respective particles becomes a certain critical point distance and generating an attractive force when the distance becomes a certain distance beyond the critical point distance, This curve or straight line is such that the curve or straight line itself changes with the temperature change, velocity change, or time change of each particle, thereby expressing a change in non-Newtonian viscosity and reproducing a change in fluid shape. particle type fluid simulation method was.

前記関数に基づいて、各粒子の圧力の経時的変化をも決定し、流体の形状の変化を再現する請求項１に記載の粒子型流体シミュレーション方法。 2. The particle-type fluid simulation method according to claim 1, wherein a time-dependent change in the pressure of each particle is also determined based on the function, and the change in the shape of the fluid is reproduced.

流体を、連続体ではなく、各微小領域を仮想的な粒子とし、各粒子と固体表面との関係において親水性及び撥水性を再現するために、各粒子と固体表面との力の相互作用を固体表面と粒子間の距離をパラメータとした関数に基づいて、各粒子の位置の経時的変化を決定し、流体と固体表面との親水性及び撥水性等の現象を再現する粒子型流体シミュレーション方法。A fluid, rather than a continuum, each micro area as virtual particles, in order to reproduce the hydrophilic and water-repellent relative to the respective particles and the solid surface, the interaction of forces between the particles and the solid surface Particle-based fluid simulation method that determines the change over time of the position of each particle based on a function using the distance between the solid surface and the particle as a parameter and reproduces phenomena such as hydrophilicity and water repellency between the fluid and the solid surface .

前記関数に基づいて、各粒子の圧力の経時的変化をも決定し、流体と固体表面との親水性及び撥水性等の現象を再現する請求項３に記載の粒子型流体シミュレーション方法。 4. The particle-type fluid simulation method according to claim 3 , wherein a time-dependent change in pressure of each particle is also determined based on the function, and phenomena such as hydrophilicity and water repellency between the fluid and the solid surface are reproduced.

前記関数は、各粒子と固体表面との間の距離がある臨界点距離となるとお互いに斥力を生じ、この臨界点距離を越えて一定の距離になるとお互いに引力を生じることを特徴とする曲線又は直線からなり、この曲線又は直線は、各粒子の温度変化、速度変化、又は時間変化に伴い、曲線又は直線自体が変化することにより、非ニュートン性粘度の変化を表現し、流体と固体表面との親水性及び撥水性等の現象を再現するようにした請求項３又は４に記載の粒子型流体シミュレーション方法。 The function is a curve characterized in that when the distance between each particle and the solid surface reaches a certain critical point distance, reciprocal forces are generated, and when the distance exceeds this critical point distance, a constant distance is generated between the particles and the solid surface. or made straight, the curve or straight line, the temperature change of each particle, the velocity change or with time variation, by curved or straight itself changes, representing the change in the non-Newtonian viscosity, fluid and solid surface 5. The method according to claim 3 , wherein phenomena such as hydrophilicity and water repellency are reproduced .