特開2023-97134シミュレーション装置、シミュレーション装置の制御方法、フィルム製造方法、塗工方法、および制御プログラム
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】公開特許公報(A)
(11)【公開番号】P2023097134
(43)【公開日】2023-07-07
(54)【発明の名称】シミュレーション装置、シミュレーション装置の制御方法、フィルム製造方法、塗工方法、および制御プログラム
(51)【国際特許分類】
G06F 30/20 20200101AFI20230630BHJP
G01N 13/02 20060101ALI20230630BHJP
B05D 3/00 20060101ALI20230630BHJP
G16Z 99/00 20190101ALI20230630BHJP
【ＦＩ】
G06F30/20
G01N13/02
B05D3/00 D
G16Z99/00
【審査請求】未請求
【請求項の数】10
【出願形態】ＯＬ
(21)【出願番号】P 2021213307
(22)【出願日】2021-12-27
(71)【出願人】
【識別番号】000002093
【氏名又は名称】住友化学株式会社
(74)【代理人】
【識別番号】110000338
【氏名又は名称】弁理士法人 ＨＡＲＡＫＥＮＺＯ ＷＯＲＬＤ ＰＡＴＥＮＴ ＆ ＴＲＡＤＥＭＡＲＫ
(72)【発明者】
【氏名】島田 直樹
(72)【発明者】
【氏名】内橋 祐介
(72)【発明者】
【氏名】八重樫 優太
(72)【発明者】
【氏名】松尾 美弥
【テーマコード（参考）】
4D075
5B146
5L049
【Ｆターム（参考）】
4D075BB91X
4D075CA35
4D075CA36
4D075CA48
4D075DC16
4D075DC22
4D075DC30
4D075EA05
5B146DJ11
5L049AA04
(57)【要約】
【課題】ＣＬＳＶＯＦ法よりも計算負荷が低く、かつ、Ｓ－ＣＬＳＶＯＦよりも高い計算精度を実現する。
【解決手段】表面張力を伴う流動計算を行うシミュレーション装置（１）は、（ｉ）空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出部（３１）と、（ｉｉ）レベルセット値算出部（３１）によって算出されたレベルセット値を用いて、注目する計算セルに含まれる界面に生じる表面張力を算出する表面張力算出部（３２）と、を備える。レベルセット値算出部（３１）は、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した計算セルから多次元的な分布情報をパラメータとして含む関数を用いてレベルセット値を算出する。
【選択図】図１
【解決手段】表面張力を伴う流動計算を行うシミュレーション装置（１）は、（ｉ）空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出部（３１）と、（ｉｉ）レベルセット値算出部（３１）によって算出されたレベルセット値を用いて、注目する計算セルに含まれる界面に生じる表面張力を算出する表面張力算出部（３２）と、を備える。レベルセット値算出部（３１）は、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した計算セルから多次元的な分布情報をパラメータとして含む関数を用いてレベルセット値を算出する。
【選択図】図１
【特許請求の範囲】
【請求項1】
表面張力を伴う流動計算を行うシミュレーション装置であって、
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出部と、
前記レベルセット値算出部によって算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出部と、を備え、
前記レベルセット値算出部は、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
シミュレーション装置。
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出部と、
前記レベルセット値算出部によって算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出部と、を備え、
前記レベルセット値算出部は、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
シミュレーション装置。
【請求項2】
前記レベルセット値算出部は、前記関数として、以下の式を用いて、前記レベルセット値（Φ）を算出する、
ここで、λは次式で与えられ、
γは次式で与えられ、
βは次式で与えられ、
ここで、βｍａｘ、βｍｉｎは予め設定された所定値、ｎｘ、ｎｙはそれぞれ界面法線のｘ方向およびｙ方向である、
請求項１に記載のシミュレーション装置。
【数1】
【数2】
【数3】
【数4】
請求項１に記載のシミュレーション装置。
【請求項3】
前記表面張力算出部は、
前記レベルセット値を用いて曲率を計算し、
前記曲率を用いて前記表面張力を算出する、
請求項１または２に記載のシミュレーション装置。
前記レベルセット値を用いて曲率を計算し、
前記曲率を用いて前記表面張力を算出する、
請求項１または２に記載のシミュレーション装置。
【請求項4】
前記表面張力算出部によって算出された前記表面張力に基づき、速度を算出する速度算出部をさらに備える、
請求項１から３のいずれか１項に記載のシミュレーション装置。
請求項１から３のいずれか１項に記載のシミュレーション装置。
【請求項5】
ＶＯＦ関数α（ｘ）を更新する更新部をさらに備え、
前記更新部は、
直前の時刻のＶＯＦ関数α（ｘ）を用いて、次式を設定し、
ここで、ｘｅｄｇｅはセル接触面の位置、Δｘはセルサイズであり、
γは次式で与えられ、
βは次式で与えられ、
ここで、βｍａｘ、βｍｉｎは予め設定された所定値、ｎｘ、ｎｙはそれぞれ界面法線のｘ方向およびｙ方向であり、
ｘＤは次式で与えられ、
ηは次式で与えられ、
前記速度算出部によって算出された前記速度（Ｖ）に基づき、前記設定したα（ｘ）のｘにｘ＋ＶΔｔを代入して、ＶＯＦ関数α（ｘ）を更新する、
請求項４に記載のシミュレーション装置。
前記更新部は、
直前の時刻のＶＯＦ関数α（ｘ）を用いて、次式を設定し、
【数5】
γは次式で与えられ、
【数6】
【数7】
ｘＤは次式で与えられ、
【数8】
【数9】
請求項４に記載のシミュレーション装置。
【請求項6】
前記表面張力算出部は、注目する計算セルに含まれる界面に生じる前記表面張力（Ｆσ）を、表面張力係数σ、界面曲率κ、および法線ベクトルＮを用いて、以下の式により計算する、
ここで、ｘ方向のセルサイズがΔｘであり、ｙ方向のセルサイズがΔｙである場合に、法線ベクトルＮのｘ成分およびｙ成分はそれぞれ、
として与えられ、
κは、
として与えられ、
ｎは、単位法線ベクトルであり、次式で与えられる、
請求項１から５のいずれか１項に記載のシミュレーション装置。
【数10】
【数11】
【数12】
κは、
【数13】
【数14】
ｎは、単位法線ベクトルであり、次式で与えられる、
【数15】
【請求項7】
表面張力を伴う流動計算を行うシミュレーション装置の制御方法であって、
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、
前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、
前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
シミュレーション装置の制御方法。
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、
前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、
前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
シミュレーション装置の制御方法。
【請求項8】
表面張力を伴う流動計算を行うシミュレーションを伴うフィルム製造方法であって、
フィルム材料の表面張力を計算するために、
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、
前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、
前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
フィルム製造方法。
フィルム材料の表面張力を計算するために、
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、
前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、
前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
フィルム製造方法。
【請求項9】
表面張力を伴う流動計算を行うシミュレーションを伴う塗工液の塗工方法であって、
前記塗工液の表面張力を計算するために、
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、
前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、
前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
塗工方法。
前記塗工液の表面張力を計算するために、
空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、
前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、
前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する、
塗工方法。
【請求項10】
請求項１に記載のシミュレーション装置としてコンピュータを機能させるための制御プログラムであって、前記レベルセット値算出部および前記表面張力算出部としてコンピュータを機能させるための制御プログラム。
【発明の詳細な説明】
【技術分野】
【0001】
本発明の一態様は、表面張力を伴う流動計算を行うシミュレーション装置等に関する。
【背景技術】
【0002】
表面張力は、薄い塗工、細い流路になるほど支配的な力であるため、近年、プリンテッドエレクトロニクス、バイオテクノロジー分野を中心としてその重要性を増している。そこで、数値計算による表面張力の予測が精力的に行われている。
【0003】
米国ロスアラモス研究センターBlackbill et. al.(1992)によって開発されたＣＳＦ（Continuum Surface Force）モデルが最もよく使用されるが、計算上の誤差が大きくなることが知られている。
【0004】
次いで、米国カルフォルニア大のOsher教授によってＬＳ（Level Set）法が開発された。その後、Sussmann教授によって、ＬＳ法とＶＯＦ（Volume of Fluid）法とを組み合わせた手法であるＣＬＳＶＯＦ（Coupling between VOF and LS）が開発された。ただし、ＣＬＳＶＯＦ法は計算負荷が大きい。
【0005】
この点を踏まえ、アイルランドのAlbadawi et al.(2013)によって、より簡便な手法である、Ｓ－ＣＬＳＶＯＦ（Simple-CLSVOF）法が開発された（非特許文献１を参照）。
【先行技術文献】
【非特許文献】
【0006】
【非特許文献1】Albadawi, A., D. B. Donoghue, A. J. Robinson, D. B. Murray, and Y. M. C. Delaure, “Influence of surface tension implementation in Volume of Fluid and coupled Volume of Fluid with Level Set methods for bubble growth and detachment,”International Journal of Multiphase Flow, 53, 11-28 (2013)
【発明の概要】
【発明が解決しようとする課題】
【0007】
Ｓ－ＣＬＳＶＯＦ法によれば、ＣＬＳＶＯＦ法に比べて計算負荷を低減できる。但し、Ｓ－ＣＬＳＶＯＦ法は、ＣＬＳＶＯＦ法に比べて、計算精度において劣っている。
【0008】
本発明の一態様は、上述の問題点に鑑みてなされたものであり、ＣＬＳＶＯＦ法よりも計算負荷が低く、かつ、Ｓ－ＣＬＳＶＯＦ法よりも高い計算精度を実現することを目的としている。
【課題を解決するための手段】
【0009】
上記の課題を解決するために、本発明の一態様に係るシミュレーション装置は、表面張力を伴う流動計算を行うシミュレーション装置であって、空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出部と、前記レベルセット値算出部によって算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出部と、を備え、前記レベルセット値算出部は、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する。
【0010】
また、本発明の一態様に係るシミュレーション装置の制御方法は、表面張力を伴う流動計算を行うシミュレーション装置の制御方法であって、空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する。
【0011】
また、本発明の一態様に係るフィルム製造方法は、表面張力を伴う流動計算を行うシミュレーションを伴うフィルム製造方法であって、フィルム材料の表面張力を計算するために、空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する。
【0012】
また、本発明の一態様に係る塗工方法は、表面張力を伴う流動計算を行うシミュレーションを伴う塗工液の塗工方法であって、前記塗工液の表面張力を計算するために、空間に設定した計算セルごとに決定された、現在の時刻の気相のセル体積占有率を用いて、注目する計算セル周りのレベルセット値を算出するレベルセット値算出ステップと、前記レベルセット値算出ステップにおいて算出された前記レベルセット値を用いて、注目する前記計算セルに含まれる界面に生じる表面張力を算出する表面張力算出ステップと、を含み、前記レベルセット値算出ステップでは、現在の時刻の気相のセル体積占有率として、ＶＯＦ（Volume of Fluid）法によって得られたセル体積占有率αＶＯＦを用いて、隣接した前記計算セルから多次元的な分布情報をパラメータとして含む関数を用いて前記レベルセット値を算出する。
【発明の効果】
【0013】
本発明の一態様によれば、ＣＬＳＶＯＦ法よりも計算負荷が低く、かつ、Ｓ－ＣＬＳＶＯＦ法よりも高い計算精度を実現できる。
【図面の簡単な説明】
【0014】
【図1】本発明の一実施形態に係るシミュレーション装置の要部の構成を示すブロック図である。
【図2】セルにおける体積率分布について説明する図である。
【図3】各セルにおけるαおよびＮの割り当てについて説明する図である。
【図4】上記シミュレーション装置において、本件手法を実行する処理の流れを例示する図である。
【図5】気泡の場合における、各手法において得られた最大速度誤差および平均速度誤差を示す図である。
【図6】液滴の場合における、各方法において得られた最大速度誤差および平均速度誤差を示す図である。
【図7】各手法によって得られた、最終時刻における気泡の速度場を例示する図である。
【図8】ＣＬＳＶＯＦ法および本件手法のそれぞれによって得られたシミュレーション結果を比較する図である。
【図9】本件手法によって得られた、各時刻における速度場を示す図である。
【発明を実施するための形態】
【0015】
以下、本発明の一実施形態に係るシミュレーション装置１について説明する。説明の便宜上、公知技術と同様の事項については、言及を適宜省略する。なお、本明細書において以下に述べる各数値は、単なる一例であることに留意されたい。また、本明細書では、数式オブジェクトの箇所を除いては、ベクトルについても、スカラーと同様の標準スタイルの文字にて表記されていることに留意されたい。
【0016】
（気液界面における動き予測についての概説）
後述する通り、本願の発明者ら（以下、単に「発明者ら」と称する）は、新たなシミュレーション手法（ＴＨＡＩＮＣ－ＬＳ法）を開発した。本明細書では、当該新たなシミュレーション手法を、本件手法とも称する。図１に示されているシミュレーション装置１の詳細な説明に先立ち、本件手法の前提となる技術的事項について簡単に述べる。
後述する通り、本願の発明者ら（以下、単に「発明者ら」と称する）は、新たなシミュレーション手法（ＴＨＡＩＮＣ－ＬＳ法）を開発した。本明細書では、当該新たなシミュレーション手法を、本件手法とも称する。図１に示されているシミュレーション装置１の詳細な説明に先立ち、本件手法の前提となる技術的事項について簡単に述べる。
【0017】
一般的に、界面における輸送に関しては、空間に設定された計算セル（計算格子）における体積率αを使用するＶＯＦ法が用いられる。体積率は、体積占有率とも称される。以下の説明では、計算セルを、単に「セル」とも略記する。
【0018】
以下、図２を参照してＶＯＦ法について説明する。図２は、セルにおける体積率分布について説明する図である。図２では、２次元空間（例：ｘｙ座標空間）におけるシミュレーションの例が示されている。但し、当業者であれば明らかである通り、以下の説明は、３次元空間（例：ｘｙｚ座標空間）におけるシミュレーションについても同様に当てはまる。
【0019】
図２の符号２０００Ａに示されている通り、気相（気体）と液相（液体）との間には、界面が存在している。図２の符号２０００Ｂには、符号２０００Ａに対応するセルにおける体積率分布を示す図である。以下の説明では、複数のセルのうち、注目する１つのセルを、注目セルと称する。
【0020】
注目セルが完全に液相（液体）で占められている場合、当該注目セルでは、α＝１である。他方、注目セルが完全に気相（気体）で占められている場合、当該注目セルでは、α＝０である。従って、注目セルにおけるαの値が０より大きく、かつ、１未満であることは、当該注目セル内に界面が存在している（言い換えれば、当該注目セルが界面を含んでいる）ことを示している。
【0021】
一例として、流体（より具体的には、非圧縮性流体）の質量保存式および運動量保存式はそれぞれ、以下の式（１）および式（２）、
…（１）
…（２）
によって与えられる。Ｖは流体の速度であり、ρは当該流体の密度であり、Ｐは当該流体の圧力であり、μは当該流体の粘度であり、ｇは重力加速度であり、Ｆσは界面に作用する表面表力である。
【数1】
【数2】
によって与えられる。Ｖは流体の速度であり、ρは当該流体の密度であり、Ｐは当該流体の圧力であり、μは当該流体の粘度であり、ｇは重力加速度であり、Ｆσは界面に作用する表面表力である。
【0022】
そして、ρおよびμはそれぞれ、以下の式（３）および式（４）、
…（３）
…（４）
によって与えられる。αＬ、ρＬ、およびμＬはそれぞれ、液相における体積率、密度、および粘度を表す。ρＬおよびμＬは、所与の物性値として与えられてよい。αＧ、ρＧ、およびμＧはそれぞれ、気相における体積率、密度、および粘度を表す。ρＧおよびμＧは、所与の物性値として与えられてよい。
【数3】
【数4】
によって与えられる。αＬ、ρＬ、およびμＬはそれぞれ、液相における体積率、密度、および粘度を表す。ρＬおよびμＬは、所与の物性値として与えられてよい。αＧ、ρＧ、およびμＧはそれぞれ、気相における体積率、密度、および粘度を表す。ρＧおよびμＧは、所与の物性値として与えられてよい。
【0023】
αは、所定の輸送方程式に従って変化する。一例として、αの輸送方程式は、以下の式（５）、
…（５）
によって与えられる。従って、式（５）を解くことにより、界面の分布および動きを予測できる。本実施形態では、シミュレーション装置１は、ＴＨＡＩＮＣ（Tangent Hyperbola with Adaptive slope for Interface Capture）法によって式（５）を解く。
【数5】
によって与えられる。従って、式（５）を解くことにより、界面の分布および動きを予測できる。本実施形態では、シミュレーション装置１は、ＴＨＡＩＮＣ（Tangent Hyperbola with Adaptive slope for Interface Capture）法によって式（５）を解く。
【0024】
図２の符号２０００Ｃには、図２の符号２０００Ｂにおける線Ｌ１に沿ったαの分布（すなわち、あるｙについての、ｘ方向におけるαの分布）を示すグラフが示されている。具体的には、符号２０００Ｃには、本件手法によって導出されたαの分布を示すグラフが示されている。当該グラフの各点は、後述のαＶＯＦに対応する。当該グラフから明らかである通り、本件手法によれば、界面の位置を正確に表現することができる。それゆえ、本件手法によれば、界面の動きを正確に予測することもできる。
【0025】
（シミュレーション装置１の構成）
続いて、図１を参照し、シミュレーション装置１の構成について述べる。図１は、シミュレーション装置１の要部の構成を示すブロック図である。以下に述べる通り、シミュレーション装置１は、本件手法によって表面張力を伴う流動計算を行うことができる。シミュレーション装置１は、制御部１０および記憶部２０を備えている。また、シミュレーション装置１には、入力装置２と表示装置３とが接続されている。
続いて、図１を参照し、シミュレーション装置１の構成について述べる。図１は、シミュレーション装置１の要部の構成を示すブロック図である。以下に述べる通り、シミュレーション装置１は、本件手法によって表面張力を伴う流動計算を行うことができる。シミュレーション装置１は、制御部１０および記憶部２０を備えている。また、シミュレーション装置１には、入力装置２と表示装置３とが接続されている。
【0026】
入力装置２は、ユーザにシミュレーション装置１を操作させるための装置である。入力装置２は、シミュレーション装置１を操作するための入力操作を受け付け、かつ、当該入力操作の内容を示す操作データをシミュレーション装置１に送信する。
【0027】
表示装置３は、シミュレーション装置１が出力する表示用データを用いて画像を表示する装置である。表示装置３には、ユーザにシミュレーション装置１を操作させるための操作画面、および、シミュレーション結果等が表示される。
【0028】
制御部１０は、シミュレーション装置１の各部を統括的に制御する。記憶部２０は、シミュレーション装置１が動作するために必要なデータおよびプログラム等を格納している。また、シミュレーションの過程において算出される途中結果、および、シミュレーションの最終結果なども、記憶部２０に格納される。
【0029】
制御部１０は、シミュレーション算出部３０、シミュレーション制御部４０、および表示制御部５０を備えている。シミュレーション算出部３０は、レベルセット値算出部３１、表面張力算出部３２、速度算出部３３、および更新部３４を備えている。シミュレーション算出部３０の各部についての説明は、後述する。
【0030】
シミュレーション制御部４０は、シミュレーション算出部３０におけるシミュレーション（より詳細には、流体シミュレーション）を制御する。一例として、シミュレーション制御部４０は、シミュレーション算出部３０に各時刻および各セルについての計算を行わせることにより、シミュレーション算出部３０にシミュレーションを実行させる。
【0031】
具体的には、シミュレーション制御部４０は、現在時刻（現在の時刻）ｔにおける各セルのαを、シミュレーション算出部３０に算出させる。以下では、現在時刻を、単に「時刻」と略記する。そして、時刻ｔにおける全てのセルにおけるαの算出が完了すると、シミュレーション制御部４０は、時刻をｔからｔ＋Δｔへと更新する。Δｔは、所定の時間刻み幅である。次いで、シミュレーション制御部４０は、時刻ｔ＋Δにおける各セルのαを、シミュレーション算出部３０に算出させる。シミュレーション制御部４０は、以上の処理を、所定の最終時刻（シミュレーション終了時刻）まで繰り返す。
【0032】
シミュレーション制御部４０は、シミュレーション算出部３０におけるシミュレーションの終了後、シミュレーション算出部３０からシミュレーション結果を取得する。一例として、シミュレーション制御部４０は、シミュレーション結果を表示制御部５０に送信してよい。あるいは、シミュレーション制御部４０は、シミュレーション結果を記憶部２０に格納してもよい。
【0033】
表示制御部５０は、シミュレーション制御部４０による上述の処理によって得られたシミュレーション結果を、表示装置３に表示させる。また、表示制御部５０は、表示させるデータを記憶部２０から読み出し、当該データに対応する画面を表示装置３に表示させてもよい。一例として、表示制御部５０は、シミュレーション装置１を操作するための操作画面データを記憶部２０から読み出し、操作画面を表示装置３に表示させる。
【0034】
（レベルセット値算出部３１）
レベルセット値算出部３１は、時刻ｔ（現在時刻）における気相のセル体積率（セル体積占有率）を用いて、注目セル周りのレベルセット値を算出する。具体的には、レベルセット値算出部３１は、現在時刻の気相のセル体積占有率として、ＶＯＦ法によって得られたセル体積占有率αＶＯＦを用いて、隣接した計算セルから多次元的な分布情報をパラメータとして含む関数を使ってレベルセット値を算出する。
レベルセット値算出部３１は、時刻ｔ（現在時刻）における気相のセル体積率（セル体積占有率）を用いて、注目セル周りのレベルセット値を算出する。具体的には、レベルセット値算出部３１は、現在時刻の気相のセル体積占有率として、ＶＯＦ法によって得られたセル体積占有率αＶＯＦを用いて、隣接した計算セルから多次元的な分布情報をパラメータとして含む関数を使ってレベルセット値を算出する。
【0035】
多次元的な分布情報をパラメータとして含む関数の例としては、（ｉ）多次元的な分布情報をパラメータとして含む高次関数、および、（ｉｉ）多次元的な分布情報をパラメータとして含む特殊関数などが挙げられる。そこで、以下では、一例として、レベルセット値算出部３１が、特殊関数を用いてレベルセット値を算出する場合を例示する。
【0036】
多次元的な分布情報をパラメータとして含む特殊関数の例としては、以下に述べる式（６）が挙げられる。なお、本発明の一態様における特殊関数とは、線形でない関数（初等関数など）を指す。特殊関数の例としては、ガンマ関数、エアリー関数、ベッセル関数、ゼータ関数、楕円関数、ルジャンドル関数、誤差関数、超幾何関数、および、直交多項式によって表される関数などが挙げられる。
【0037】
本実施形態では、レベルセット値算出部３１は、多次元的な分布情報をパラメータとして含む関数（より具体的には、特殊関数）として、以下の式（６）を用いて、レベルセット値（Φ）を算出する、
…（６）
ここで、λは以下の式（７）によって与えられ、
…（７）
γは以下の式（８）によって与えられ、
…（８）
βは以下の式（９）によって与えられ、
…（９）
ここで、βｍａｘ、βｍｉｎは予め設定された所定値、ｎｘ、ｎｙはそれぞれ界面法線のｘ方向およびｙ方向である。なお、式（６）および（７）の導出過程については後述する。
【数6】
ここで、λは以下の式（７）によって与えられ、
【数7】
γは以下の式（８）によって与えられ、
【数8】
βは以下の式（９）によって与えられ、
【数9】
ここで、βｍａｘ、βｍｉｎは予め設定された所定値、ｎｘ、ｎｙはそれぞれ界面法線のｘ方向およびｙ方向である。なお、式（６）および（７）の導出過程については後述する。
【0038】
（表面張力算出部３２）
表面張力算出部３２は、レベルセット値算出部３１によって算出されたレベルセット値（例：上述のΦ）を用いて、注目セルに含まれる界面に生じる表面張力（Ｆσ）を算出する。
表面張力算出部３２は、レベルセット値算出部３１によって算出されたレベルセット値（例：上述のΦ）を用いて、注目セルに含まれる界面に生じる表面張力（Ｆσ）を算出する。
【0039】
一例として、表面張力算出部３２は、レベルセット値を用いて曲率を計算し、当該曲率を用いて表面張力を算出してよい。
【0040】
本実施形態では、表面張力算出部３２は、以下の式（１０）、
…（１０）
によって、Ｆσを算出する、
ここで、σは表面張力係数であり、κは界面曲率であり、Ｎは法線ベクトルであり、
Ｎのｘ成分およびｙ成分はそれぞれ、以下の式（１１）～（１２）、
…（１１）
…（１２）
によって与えられ、
ここで、Δｘはｘ方向のセルサイズであり、Δｙはｙ方向のセルサイズであり、
κは、以下の式（１３）～（１４）、
…（１３）
…（１４）
によって与えられ、
ｎは単位法線ベクトルであり、以下の式（１５）によって与えられる。
…（１５）
【数10】
によって、Ｆσを算出する、
ここで、σは表面張力係数であり、κは界面曲率であり、Ｎは法線ベクトルであり、
Ｎのｘ成分およびｙ成分はそれぞれ、以下の式（１１）～（１２）、
【数11】
【数12】
によって与えられ、
ここで、Δｘはｘ方向のセルサイズであり、Δｙはｙ方向のセルサイズであり、
κは、以下の式（１３）～（１４）、
【数13】
【数14】
によって与えられ、
ｎは単位法線ベクトルであり、以下の式（１５）によって与えられる。
【数15】
【0041】
以上の通り、表面張力算出部３２は、Φを用いることにより、σ、κ、Ｎ、Δｘ、およびΔｙを用いて、Ｆσを算出できる。σは、所与の物性値として与えられてよい。ΔｘおよびΔｙは、シミュレーション装置１のユーザによって設定されてよい。なお、式（１０）は、ＣＳＦモデルに基づくＦσの算出式の一例である。
【0042】
式（１１）～（１２）に関し、図３を参照して説明する。図３は、各セルにおけるαおよびＮの割り当てについて説明する図である。ｉはｘ方向のセル位置を示す添字であり、ｊはｙ方向のセル位置を表す添字である。本明細書では、あるセル（ｉ，ｊ）におけるαを、αｉ，ｊと表記する。
【0043】
図３に示される通り、αは、各セルの中心に対応付けられている。他方、式（１１）および（１２）によれば、Ｎは２つの隣り合うセル間の境界に対応付けられる。式（１１）におけるＮｘ，ｉ＋１／２，ｊは、セル（ｉ，ｊ）とセル（ｉ＋１，ｊ）との境界における、Ｎのｘ成分である。また、式（１２）におけるＮｙ，ｉ，ｊ＋１／２は、セル（ｉ，ｊ）とセル（ｉ，ｊ＋１）との境界における、Ｎのｙ成分である。
【0044】
式（１４）は、ＣＳＦモデルに基づくκの算出式の一例である。式（１４）によれば、各セルの中心におけるκを算出できる。式（１５）は、単位法線ベクトルの算出式の一例である。式（１３）は、κに関する公知の補間式の一例である。式（１３）によれば、２つの隣り合うセル間の境界におけるκを算出できる。
【0045】
（速度算出部３３）
速度算出部３３は、表面張力算出部３２によって算出された表面張力に基づき、速度Ｖを算出する。一例として、速度算出部３３は、公知の算出式を用いて、Ｆσに基づき、Ｖを算出してよい。
速度算出部３３は、表面張力算出部３２によって算出された表面張力に基づき、速度Ｖを算出する。一例として、速度算出部３３は、公知の算出式を用いて、Ｆσに基づき、Ｖを算出してよい。
【0046】
なお、シミュレーション算出部３０は、公知の算出式を用いて、Ｐを算出してもよい。本実施形態では、説明の便宜上、速度算出部３３にＰを算出する機能が付与されている場合を例示する。
【0047】
（更新部３４）
更新部３４は、ＶＯＦ関数α（ｘ）を更新する。まず、更新部３４は、直前の時刻（前タイムステップにおける時刻，すなわちｔ－Δｔ）におけるα（ｘ）を用いて、以下の式（１６）、
…（１６）
を設定する。
更新部３４は、ＶＯＦ関数α（ｘ）を更新する。まず、更新部３４は、直前の時刻（前タイムステップにおける時刻，すなわちｔ－Δｔ）におけるα（ｘ）を用いて、以下の式（１６）、
【数16】
を設定する。
【0048】
式（１６）において、ｘｅｄｇｅはセル接触面の位置であり、Δｘはセルサイズ（より具体的には、ｘ方向のセルサイズ）である。γは、上述の式（８）によって与えられる。βは、上述の式（９）によって与えられる。
ここで、ｘＤは、以下の式（１７）によって与えられ、
…（１７）
ηは、以下の式（１８）によって与えられる。
…（１８）
ここで、ｘＤは、以下の式（１７）によって与えられ、
【数17】
ηは、以下の式（１８）によって与えられる。
【数18】
【0049】
そして、更新部３４は、速度算出部３３によって算出されたＶに基づき、上述の通り設定したα（ｘ）のｘにｘ＋ＶΔｔを代入して、現在時刻のα（ｘ）を算出する（すなわち、α（ｘ）を更新する）。このように、更新部３４は、Δｔの経過に伴ってｘをＶΔｔだけ移動させることにより、α（ｘ）を更新する。
【0050】
（本件手法についての補足説明）
発明者らは、上述の式（１６）に基づき、以下の式（１９）、
…（１９）
を新たに導出した。式（１９）は、式（１６）の原点を各セルの中心に設定し（すなわち、ｘｅｄｇｅを０に設定し）、かつ、ｘＤをΦに置き換えることによって得られる。
発明者らは、上述の式（１６）に基づき、以下の式（１９）、
【数19】
を新たに導出した。式（１９）は、式（１６）の原点を各セルの中心に設定し（すなわち、ｘｅｄｇｅを０に設定し）、かつ、ｘＤをΦに置き換えることによって得られる。
【0051】
さらに、発明者らは、α（ｘ）について、以下の式（２０）、
…（２０）
によって示される条件を与えるという新たな着想に至った。
【数20】
によって示される条件を与えるという新たな着想に至った。
【0052】
発明者らは、式（１９）を式（２０）に代入して式変形することによって、上述の式（６）～（７）を導出した。このように、式（６）～（７）は、発明者らによって新たに提案された式である。
【0053】
上述の通り、本件手法は、ＴＨＡＩＮＣ法に基づいている。加えて、本件手法では、α（ｘ）を積分したαＶＯＦを用いて、レベルセット値（ＬＳ値）Φが定義されている（上述の式（６）～（７）および（２０）を参照）。このことから、発明者らは、本件手法を、ＴＨＡＩＮＣ－ＬＳ法とネーミングした。
【0054】
（シミュレーション装置１における処理の流れ）
図４は、シミュレーション装置１において、本件手法を実行する処理の流れを例示する図である。なお、図４の例では、シミュレーション装置１における処理の開始に先立ち、ρＬ、ρＧ、μＬ、μＧ、およびσが、所与の物性値として与えられているものとする。
図４は、シミュレーション装置１において、本件手法を実行する処理の流れを例示する図である。なお、図４の例では、シミュレーション装置１における処理の開始に先立ち、ρＬ、ρＧ、μＬ、μＧ、およびσが、所与の物性値として与えられているものとする。
【0055】
一例として、シミュレーション制御部４０は、ユーザによる所定の入力操作を受け付けて、シミュレーション算出部３０にシミュレーションを開始させる。まず、シミュレーション制御部４０は、シミュレーション算出部３０にシミュレーション条件を初期化させる（Ｓ１）。
【0056】
初期化が完了した後、シミュレーション算出部３０は、αの分布からρおよびμを算出する（Ｓ２）。具体的には、シミュレーション算出部３０は、上述の式（３）および（４）に従って、αに応じたρおよびμを算出する。
【0057】
次いで、シミュレーション算出部３０は、αを用いてΦを算出する（Ｓ３）。具体的には、上述の通り、レベルセット値算出部３１は、αＶＯＦを用いてΦを算出する。Ｓ３は、レベルセット値算出ステップと称されてもよい。
【0058】
次いで、シミュレーション算出部３０は、Ｆσを算出する（Ｓ４）。具体的には、上述の通り、表面張力算出部３２は、Ｓ３において算出されたΦを用いて、Ｆσを算出する。Ｓ４は、表面張力算出ステップと称されてもよい。
【0059】
次いで、シミュレーション算出部３０は、ＶおよびＰを算出する（Ｓ５）。具体的には、上述の通り、速度算出部３３は、Ｓ４において算出されたＦσに基づき、Ｖを算出する。そして、速度算出部３３は、Ｐをさらに算出する。
【0060】
次いで、シミュレーション算出部３０は、αを更新する（Ｓ６）。具体的には、上述の通り、更新部３４は、直前時刻におけるα（ｘ）を用いて、速度算出部３３によって算出されたＶに基づき、α（ｘ）を更新する。
【0061】
次いで、シミュレーション制御部４０は、時刻ｔが最終時刻（シミュレーション終了時刻）に達したか否かを判定する（Ｓ７）。時刻ｔが最終時刻に達した場合（Ｓ７においてＹＥＳ）、シミュレーション制御部４０は、シミュレーション算出部３０にシミュレーションを終了させる。
【0062】
他方、時刻ｔが最終時刻に達していない場合（Ｓ７においてＮＯ）、シミュレーション制御部４０は、時刻ｔを更新する（Ｓ８）。具体的には、シミュレーション制御部４０は、ｔをｔ＋Δｔに更新する。そして、Ｓ２に戻る。このように、初期時刻から最終時刻に至るまで、一連の処理が繰り返される。図４の例では、シミュレーションの終了後、シミュレーション装置１は、初期時刻から最終時刻までの各時刻におけるα、Ｖ、およびＰを、シミュレーション結果として出力してよい。
【0063】
（検証結果）
発明者らは、本件手法（ＴＨＡＩＮＣ－ＬＳ法）の有効性を検証するために、複数の数値シミュレーションを行った。以下、一部の検証結果について述べる。
発明者らは、本件手法（ＴＨＡＩＮＣ－ＬＳ法）の有効性を検証するために、複数の数値シミュレーションを行った。以下、一部の検証結果について述べる。
【0064】
（ＶＯＦ法およびＳ－ＣＬＳＶＯＦとの比較）
まず、発明者らは、気泡および液滴のそれぞれについて、ＶＯＦ法、Ｓ－ＣＬＳＶＯＦ法、および本件手法のそれぞれによってシミュレーションを行うことにより、各方法を比較した。具体的には、発明者らは、各手法によって得られた最大速度誤差｜Ｖｍａｘ｜および平均速度誤差｜Ｖａｖｅ｜を比較した。
まず、発明者らは、気泡および液滴のそれぞれについて、ＶＯＦ法、Ｓ－ＣＬＳＶＯＦ法、および本件手法のそれぞれによってシミュレーションを行うことにより、各方法を比較した。具体的には、発明者らは、各手法によって得られた最大速度誤差｜Ｖｍａｘ｜および平均速度誤差｜Ｖａｖｅ｜を比較した。
【0065】
なお、最大速度誤差および平均速度誤差はそれぞれ、以下の式（２１）および（２２）、
…（２１）
…（２２）
によって与えられる。ここで、Ｖｉ，ｊは、セル（ｉ，ｊ）におけるＶを表す。
【数21】
【数22】
によって与えられる。ここで、Ｖｉ，ｊは、セル（ｉ，ｊ）におけるＶを表す。
【0066】
発明者らは、Ｃｏａｒｓｅ（粗）およびＦｉｎｅ（密）という２つの異なるセル条件を設定し、各セルサイズにおける最大速度誤差および平均速度誤差を、各手法によって得た。Ｃｏａｒｓｅにおいて、セル数は５０×５０＝２５００であり、セルサイズは、Δｘ＝Δｙ＝０．００１ｍである。他方、Ｆｉｎｅにおいて、セル数は１００×１００＝１００００であり、セルサイズは、Δｘ＝Δｙ＝０．０００５ｍである。
【0067】
【0068】
図５および図６に示される通り、Ｓ－ＣＬＳＶＯＦ法では、多くのケースにおいて、ＶＯＦ法よりも小さい最大速度誤差および平均速度誤差が得られた。そして、図５および図６に示される通り、本件手法によれば、気泡および液滴のいずれについても、かつ、ＣｏａｒｓｅおよびＦｉｎｅのいずれについても、Ｓ－ＣＬＳＶＯＦ法よりも小さい最大速度誤差および平均速度誤差が得られた。このように、本件手法によって、Ｓ－ＣＬＳＶＯＦ法よりも高い計算精度を実現できることが裏付けられた。
【0069】
図７には、各手法によって得られた、最終時刻における気泡の速度場（より具体的には、Ｆｉｎｅにおける速度場）が示されている。図７において、符号７０００ＡはＶＯＦ法における速度場を、符号７０００ＢはＳ－ＣＬＳＶＯＦ法における速度場を、符号７０００Ｃは本件手法における速度場を、それぞれ示す。
【0070】
図７に示される通り、Ｓ－ＣＬＳＶＯＦ法によれば、ＶＯＦ法において顕著に生じていたスプリアス流（spurious current）を低減できる。そして、本件手法によれば、Ｓ－ＣＬＳＶＯＦ法よりも効果的にスプリアス流低減できることがさらに確認された。このことからも、本件手法によって、Ｓ－ＣＬＳＶＯＦ法よりも高い計算精度を実現できることが裏付けられた。
【0071】
なお、発明者らは、本件手法に従ってΦを定義することにより、Ｓ－ＣＬＳＶＯＦ法に比べてΦの精度を向上させることができると考えている。このことから、発明者らは、上述した本件手法の高い計算精度は、Φの精度向上によってもたらされたと考えている。
【0072】
（ＣＬＳＶＯＦ法との比較）
続いて、発明者らは、ＣＬＳＶＯＦ法および本件手法のそれぞれによって、ノズルからの気泡生成をシミュレートした。図８は、ＣＬＳＶＯＦ法および本件手法のそれぞれによって得られたシミュレーション結果を比較する図である。図８には、各時刻（ｔ１、ｔ２、ｔ３、およびｔ４）における、ＣＬＳＶＯＦ法による結果（左側）と本件手法による結果（右側）とが示されている。
続いて、発明者らは、ＣＬＳＶＯＦ法および本件手法のそれぞれによって、ノズルからの気泡生成をシミュレートした。図８は、ＣＬＳＶＯＦ法および本件手法のそれぞれによって得られたシミュレーション結果を比較する図である。図８には、各時刻（ｔ１、ｔ２、ｔ３、およびｔ４）における、ＣＬＳＶＯＦ法による結果（左側）と本件手法による結果（右側）とが示されている。
【0073】
図８の例では、ｔ１において１番目の気泡（最初の気泡）が生成されており、かつ、ｔ４において２番目の気泡（後続する気泡）が生成されている。図８の例は、最初の気泡が、後続する気泡の生成に不安定な影響を及ぼしていないことを表している。
【0074】
図８に示される通り、いずれの時刻においても、本件手法による結果は、ＣＬＳＶＯＦ法による結果とほぼ一致することが確認された。このことから、本件手法によれば、ＣＬＳＶＯＦ法に近い計算精度（言い換えれば、Ｓ－ＣＬＳＶＯＦ法よりも高い計算精度）を実現できることが裏付けられた。
【0075】
なお、本件手法によれば、上述の通り、比較的簡便な式によってΦを定義できる。このことから、本件手法によれば、ＣＬＳＶＯＦ法よりも低い計算負荷によってΦを定義できる。これらのことから、本件手法は、ＣＬＳＶＯＦ法よりも計算負荷が低く、かつ、Ｓ－ＣＬＳＶＯＦよりも高い計算精度を実現できる手法であると言える。
【0076】
図９には、本件手法によって得られた、上述の各時刻における速度場が示されている。具体的には、図９には、ピンチオフ時の速度場の時刻歴が示されている。図９の例は、気泡離脱直前においてガスの狭流域が急速に形成され、当該狭流域が急速な速度上昇と相互干渉する様子を示している。図９に示されている各時刻における速度場は、図８に示されている、各時刻における本件手法の結果と良好に対応している。
【0077】
（効果）
以上の通り、シミュレーション装置１によれば、本件手法を実行できる。具体的には、シミュレーション装置１では、レベルセット値算出部３１は、αＶＯＦを用いて、隣接したセルから多次元的な分布情報をパラメータとして含む関数を用いて、レベルセット値（例：Φ）を算出する。当該構成によれば、２次元および３次元上の勾配を有する界面におけるレベルセット値をより正確に算出することができ、かつ、反復計算による計算負荷が低い状態にて当該レベルセット値を算出することが可能となる。
以上の通り、シミュレーション装置１によれば、本件手法を実行できる。具体的には、シミュレーション装置１では、レベルセット値算出部３１は、αＶＯＦを用いて、隣接したセルから多次元的な分布情報をパラメータとして含む関数を用いて、レベルセット値（例：Φ）を算出する。当該構成によれば、２次元および３次元上の勾配を有する界面におけるレベルセット値をより正確に算出することができ、かつ、反復計算による計算負荷が低い状態にて当該レベルセット値を算出することが可能となる。
【0078】
従って、例えば、上述の通り、本件手法によれば、Ｓ－ＣＬＳＶＯＦ法と同程度の計算負荷にて、ＣＬＳＶＯＦ法に近い計算精度が得られる。すなわち、本件手法によれば、ＣＬＳＶＯＦ法よりも計算負荷が低く、かつ、Ｓ－ＣＬＳＶＯＦよりも高い計算精度を実現できる。
【0079】
（シミュレーション装置１の応用例）
流体には、様々な種類の力が作用する。これらの力の例としては、慣性力、重力、圧力、および粘性力に加えて、表面張力を挙げることができる。流体力学の分野における、表面張力に関する無次元数としては、例えば、キャピラリー数（粘性力と表面張力との比）、ボンド数（重力と表面張力との比，エトベス数とも称される）、およびウェーバー数（慣性力と表面張力との比）が知られている。
流体には、様々な種類の力が作用する。これらの力の例としては、慣性力、重力、圧力、および粘性力に加えて、表面張力を挙げることができる。流体力学の分野における、表面張力に関する無次元数としては、例えば、キャピラリー数（粘性力と表面張力との比）、ボンド数（重力と表面張力との比，エトベス数とも称される）、およびウェーバー数（慣性力と表面張力との比）が知られている。
【0080】
これらの無次元数が１を下回る場合には、運動方程式において、表面張力がその他の力に比べて支配的であると考えられる。すなわち、これらの無次元数が１を下回る場合には、表面張力を無視することはできないと考えられる。
【0081】
工学上の様々な分野において、表面張力を考慮した設計が必要とされる。例えば、撥水加工、ハードディスクドライブＨＤＤ内のベアリングの摩擦低減、接着、内視鏡または血管手術用のステントワイヤーの製造、コンタクトレンズの製造、および電子回路の製造のために、表面張力を考慮した設計が必要とされる。
【0082】
多くの流体では、概ね数１０ｎｍ～０．１ｍｍまたはそれより小さいスケールにおける現象が考察されている。そこで、シミュレーション装置１を用いて、当該スケールにおける現象を高精度にシミュレーションすることにより、表面張力を考慮した設計における誤差を低減することが可能となる。一例として、シミュレーション装置１によるシミュレーションは、フィルム製造方法に適用されてよい。
【0083】
従って、例えば、本発明の一態様に係るフィルム製造方法は、表面張力を伴う流動計算を行うシミュレーションを伴いうる。当該フィルム製造方法は、フィルム材料の表面張力を計算するために、上述のレベルセット値算出ステップと表面張力算出ステップとを含んでいてよい。
【0084】
別の例として、シミュレーション装置１によるシミュレーションは、塗工液の塗工方法に適用されてもよい。本発明の一態様に係る塗工方法も、表面張力を伴う流動計算を行うシミュレーションを伴いうる。当該塗工方法は、塗工液の表面張力を計算するために、上述のレベルセット値算出ステップと表面張力算出ステップとを含んでいてよい。
【0085】
〔ソフトウェアによる実現例〕
シミュレーション装置１（以下、「装置」と呼ぶ）の機能は、当該装置としてコンピュータを機能させるためのプログラムであって、当該装置の各制御ブロック（特に制御部１０に含まれる各部）としてコンピュータを機能させるためのプログラムにより実現することができる。
シミュレーション装置１（以下、「装置」と呼ぶ）の機能は、当該装置としてコンピュータを機能させるためのプログラムであって、当該装置の各制御ブロック（特に制御部１０に含まれる各部）としてコンピュータを機能させるためのプログラムにより実現することができる。
【0086】
この場合、上記装置は、上記プログラムを実行するためのハードウェアとして、少なくとも１つの制御装置（例えばプロセッサ）と少なくとも１つの記憶装置（例えばメモリ）を有するコンピュータを備えている。この制御装置と記憶装置により上記プログラムを実行することにより、上記実施形態で説明した各機能が実現される。
【0087】
上記プログラムは、一時的ではなく、コンピュータ読み取り可能な、１または複数の記録媒体に記録されていてもよい。この記録媒体は、上記装置が備えていてもよいし、備えていなくてもよい。後者の場合、上記プログラムは、有線または無線の任意の伝送媒体を介して上記装置に供給されてもよい。
【0088】
また、上記各制御ブロックの機能の一部または全部は、論理回路により実現することも可能である。例えば、上記各制御ブロックとして機能する論理回路が形成された集積回路も本発明の一態様の範疇に含まれる。この他にも、例えば量子コンピュータにより上記各制御ブロックの機能を実現することも可能である。
【0089】
また、上記実施形態で説明した各処理は、ＡＩ（Artificial Intelligence：人工知能）に実行させてもよい。この場合、ＡＩは上記制御装置で動作するものであってもよいし、他の装置（例えばエッジコンピュータまたはクラウドサーバ等）で動作するものであってもよい。
【0090】
〔付記事項〕
本発明の一態様は上述した実施形態に限定されるものではなく、請求項に示した範囲で種々の変更が可能であり、異なる実施形態にそれぞれ開示された技術的手段を適宜組み合わせて得られる実施形態についても本発明の一態様の技術的範囲に含まれる。
本発明の一態様は上述した実施形態に限定されるものではなく、請求項に示した範囲で種々の変更が可能であり、異なる実施形態にそれぞれ開示された技術的手段を適宜組み合わせて得られる実施形態についても本発明の一態様の技術的範囲に含まれる。
【符号の説明】
【0091】
１ シミュレーション装置
１０ 制御部
３０ シミュレーション算出部
３１ レベルセット値算出部
３２ 表面張力算出部
３３ 速度算出部
３４ 更新部
４０ シミュレーション制御部
１０ 制御部
３０ シミュレーション算出部
３１ レベルセット値算出部
３２ 表面張力算出部
３３ 速度算出部
３４ 更新部
４０ シミュレーション制御部
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
