(57)【要約】【課題】第１流体、第２流体、および第１流体と第２流体との間のインタフェースの動きをシミュレートするプロセス。【解決手段】第１流体および第２流体の流体速度および圧力が第１の時間ステップにおいて第１メッシュ上の第１組のノードにおいて計算される。流体速度および圧力は第１の時間ステップにおける第１流体および第２流体の動きを表す。レベルセットは第１の時間ステップにおけるインタフェースの位置を表すが、第２メッシュにおける第２組のノードにおいて計算される。第２メッシュは第１メッシュより奇数倍細かい。第１流体および第２流体の新規流体速度および圧力が第２の時間ステップにおいて第１組の点において計算される。第１流体および第２流体の新流体速度および圧力は第２の時間ステップにおける第１組のノードで計算されるが、第２の時間ステップにおける第１流体および第２流体の新規の動きを表す。【選択図】図８
2008047109.tif 000002
【特許請求の範囲】
【請求項１】
第１流体、第２流体、および前記第１流体と前記第２流体との間のインタフェースの動
きをシミュレートするプロセスで、前記第１流体および前記第２流体の流体速度および圧
力を、第１の時間ステップにおいて第１メッシュ上の第１組のノードで計算するステップ
で、これにより前記第１の時間ステップにおける前記第１流体および前記第２流体の動き
を表すステップと、前記第１の時間ステップにおける前記インタフェースの位置を表すレ
ベルセットを、前記第１メッシュより奇数倍細かい第２メッシュ上の第２組のノードで計
算するステップと、前記第１流体および前記第２流体の新しい流体速度および圧力を、第
２の時間ステップにおいて前記第１組のノードで計算するステップで、これにより前記第
２の時間ステップにおける前記第１流体および前記第２流体の新しい動きを表すステップ
と、前記第２組のノードにおける新しいレベルセットを計算するステップで、これにより
前記第２の時間ステップにおける前記インタフェースの新しい位置を表すステップとを含
むプロセス。
【請求項２】
前記第２メッシュは前記第１メッシュより３倍細かい、請求項１に記載されるプロセス
。
【請求項３】
前記流体速度および圧力の計算は有限差分方法で行なわれる、請求項１に記載されるプ
ロセス。
【請求項４】
前記第１流体はインクで前記第２流体は空気である、請求項１に記載されるプロセス。
【請求項５】
シミュレートされる前記動きは前記第１流体の一部がノズルから射出されることを含む
、請求項１に記載されるプロセス。
【請求項６】
前記第１メッシュおよび前記第２メッシュは物理的空間における四辺形メッシュである
、請求項１に記載されるプロセス。
【請求項７】
さらに、前記第１メッシュおよび前記第２メッシュを前記物理的空間から、前記第１メ
ッシュおよび前記第２メッシュが長方形の演算空間にマッピングすること、を含む、請求
項６に記載されるプロセス。
【請求項８】
さらに前記物理的空間を前記演算空間にマッピングするためにプロセスの初めに前記第
１メッシュおよび前記第２メッシュ向けに２組の変換マトリクスおよびヤコビアンを計算
し格納するステップを含む、請求項７に記載されるプロセス。
【請求項９】
前記流体の前記動きは三次元の座標システムにおいてシミュレートされる、請求項１に
記載されるプロセス。
【請求項１０】
前記座標システムは軸対称座標システムで、方位に沿った角運動はシミュレートされな
い、請求項９に記載されるプロセス。
【請求項１１】
前記流体速度は前記第１メッシュの中心において計算され、前記圧力は前記第１メッシ
ュのノードにおいて計算され、前記レベレセットは前記第２メッシュの前記ノードにおい
て計算される、請求項１に記載されるプロセス。
【請求項１２】
さらに半時間ステップにおいて前記第２メッシュのエッジ上におけるエッジ速度の計算
を含み、
前記半時間ステップは前記第１の時間ステップおよび前記第２の時間ステップの中間に
あり、前記エッジのレベルセット値は、前記第１メッシュのエッジと一致する前記第２メ
ッシュのエッジに対してはテイラー展開を用い、また前記第２メッシュ上のエッジ速度に
対しては距離加重平均を用いて計算される、請求項１に記載されるプロセス。
【請求項１３】
さらに、第１軸に沿った前記第２組のノードにおける前記レベルセットの中央差の計算
と、前記第１軸に沿った前記第２組のノードにおける前記レベルセットの前方差の計算と
、前記第１軸に沿った前記第２組のノードにおける前記レベルセットの後方差の計算と、
第１軸に沿った前記第２組のノードにおける前記レベルセットの限定的導関数の計算と、
第１軸に沿った前記第２組のノードにおける二次の限定的導関数の計算と、第２軸に沿っ
た第２組のノードにおける接線導関数項の計算、を含む、請求項１に記載されるプロセス
。
【請求項１４】
装置が請求項１に記載されるプロセスを実行するよう指示する１組の命令を含むコンピ
ュータ・プログラムをコード化した機械読み取り可能な媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明はレベルセット方法を用いて２つの流体の動きおよびこれらの流体間にインタフ
ェースをモデル化、シミュレート、および分析するシステムおよび方法に関する。
【背景技術】
【０００２】
インクジェット印刷ヘッドは印刷媒体にインクの液滴を射出することにより画像を実現
する印刷装置である。インク射出プロセスとそれに続くインクの液滴を制御することは印
刷ヘッドにより作成される任意の産物の質を確保するために不可欠である。このような制
御を達成するためには印刷および射出プロセスの正確かつ効率的なシミュレーションを行
なうことが重要である。このプロセスをシミュレートすることは少なくとも２つの流体（
すなわちインクと空気）およびこれら流体間のインタフェースのモデル化を含む。先行技
術による方法はこの挙動をモデル化するために計算流体力学、有限要素解析、有限差分解
析、およびレベルセット方法を用いている。
【０００３】
レベルセット方法は有用であるが、質量保存について問題を有することがある。質量保
存とはシミュレーション方法において閉鎖システムの場合シミュレートされるシステムの
質量を保存する能力を指す。同様に開放システムにおいて、増減された質量はシステムの
シミュレートされる質量に反映される必要がある。先行技術における方法は定期的にレベ
ルセットを再距離化するか、またはシミュレーション全体に対しより細かいメッシュを用
いることにより質量保存の問題に対応している。
【０００４】
【特許文献１】米国特許出願公開第２００４０１８１３８３号明細書
【発明の開示】
【発明が解決しようとする課題】
【０００５】
より細かいメッシュを使用することはシミュレーションの質量保存を改良する一方、シ
ミュレーションに必要な時間およびリソースが甚だしい増大につながる。例えば、二次元
シミュレーションに対しｎ×ｍメッシュで始めた場合２×ｎ×ｍの速度、ｎ×ｍの圧力、
およびｎ×ｍのレベルセットを格納する必要がある。流体の非圧縮性を実施するために、
ｎ×ｍ掛けるｎ×ｍの係数メトリクスを有する線形システムを解く必要がある。このよう
に、メッシュが２倍細かくなると、格納要求は４倍になる。さらに、非ゼロ要素を４倍含
む係数マトリクスを有する線形システムを解く必要がある。直接ソルバが用いられる場合
、４倍多い係数マトリクスは通常６４倍長い演算時間を必要とする。これはさらに必要な
演算時間を増加させることになる。
【０００６】
従って、演算リソースに不必要に負担をかけずにレベルセット方法の質量保存を改良す
るシミュレーション方法が必要とされる。
【０００７】
本発明の一目的はレベルセット方法を用いるシミュレーションにおいて、シミュレーシ
ョンに必要な演算リソースを制限しながら質量保存を改良することである。
【課題を解決するための手段】
【０００８】
本発明の一実施形態は第１流体、第２流体、および第１流体と第２流体と間のインタフ
ェースの動きをシミュレートするシステムまたは方法であることができる。第１流体およ
び第２流体の流体速度および圧力が第１の時間ステップにおいて第１メッシュ上の第１組
の点において計算される。流体速度および圧力は第１の時間ステップにおける第１流体お
よび第２流体の動きを表す。レベルセットは第１の時間ステップにおけるインタフェース
の位置を表すが、第１メッシュより奇数倍細かい第２メッシュにおける第２組の点におい
て計算される。第１流体および第２流体の新規流体速度および圧力が第２の時間ステップ
において第１組の点において計算される。新規流体速度および圧力は第２の時間ステップ
における第１流体および第２流体の新規の動きを表す。第２組の点における新規レベルセ
ットが計算される。新規レベルセットはインタフェースの第２の時間ステップにおける新
位置を表す。
【０００９】
添付図面と併せ、以下の説明およびクレームを参照することにより他の目的や成果と共
に発明のより完全な理解が明らかになり、認識されよう。
【００１０】
本発明の、第１流体、第２流体、および前記第１流体と前記第２流体との間のインタフ
ェースの動きをシミュレートするプロセスは、前記第１流体および前記第２流体の流体速
度および圧力を、第１の時間ステップにおいて第１メッシュ上の第１組のノードで計算す
るステップで、これにより前記第１の時間ステップにおける前記第１流体および前記第２
流体の動きを表すステップと、前記第１の時間ステップにおける前記インタフェースの位
置を表すレベルセットを、前記第１メッシュより奇数倍細かい第２メッシュ上の第２組の
ノードで計算するステップと、前記第１流体および前記第２流体の新しい流体速度および
圧力を、第２の時間ステップにおいて前記第１組のノードで計算するステップで、これに
より前記第２の時間ステップにおける前記第１流体および前記第２流体の新しい動きを表
すステップと、前記第２組のノードにおける新しいレベルセットを計算するステップで、
これにより前記第２の時間ステップにおける前記インタフェースの新しい位置を表すステ
ップとを含む。
【００１１】
また、前記第２メッシュは前記第１メッシュより３倍細かくてもよい。
【００１２】
また、前記流体速度および圧力の計算は有限差分方法で行なわれるものであってもよい
。
【００１３】
また、前記第１流体はインクで前記第２流体は空気であってもよい。
【００１４】
また、シミュレートされる前記動きは前記第１流体の一部がノズルから射出されること
を含むものであってもよい。
【００１５】
また、前記第１メッシュおよび前記第２メッシュは物理的空間における四辺形メッシュ
であってもよい。
【００１６】
また、前記第１メッシュおよび前記第２メッシュを前記物理的空間から、前記第１およ
び前記第２メッシュが長方形の演算空間にマッピングすることを含むものであってもよい
。
【００１７】
また、さらに前記物理的空間を前記演算空間にマッピングするためにプロセスの初めに
前記第１メッシュおよび前記第２メッシュ向けに２組の変換マトリクスおよびヤコビアン
を計算し格納するステップを含むものであってもよい。
【００１８】
また、前記流体の前記動きは三次元の座標システムにおいてシミュレートされるもので
あってもよい。
【００１９】
また、前記座標システムは軸対称座標システムで、方位に沿った角運動はシミュレート
されないものであってもよい。
【００２０】
また、前記流体速度は前記第１メッシュの中心において計算され、前記圧力は前記第１
メッシュのノードにおいて計算され、前記レベレセットは前記第２メッシュの前記ノード
において計算されるものであってもよい。
【００２１】
また、さらに半時間ステップにおいて前記第２メッシュのエッジ上におけるエッジ速度
の計算を含み、前記半時間ステップは前記第１の時間ステップおよび前記第２の時間ステ
ップの中間にあり、前記エッジのレベルセット値は、前記第１メッシュのエッジと一致す
る前記第２メッシュのエッジに対してはテイラー展開を用い、また前記第２メッシュ上の
エッジ速度に対しては距離加重平均を用いて計算されるものであってもよい。
【００２２】
また、さらに、第１軸に沿った前記第２組のノードにおける前記レベルセットの中央差
の計算と、前記第１軸に沿った前記第２組のノードにおける前記レベルセットの前方差の
計算と、前記第１軸に沿った前記第２組のノードにおける前記レベルセットの後方差の計
算と、第１軸に沿った前記第２組のノードにおける前記レベルセットの限定的導関数（a 
limiting derivative）の計算と、第１軸に沿った前記第２組のノードにおける二次の限
定的（a second-order limiting derivative）導関数の計算と、第２軸に沿った第２組の
ノードにおける接線導関数項の計算を含むものであってもよい。
【００２３】
また、本発明は、上述のプロセスを実行するよう指示する１組の命令を含むコンピュー
タ・プログラムをコード化した機械読み取り可能な媒体であってもよい。
【発明を実施するための最良の形態】
【００２４】
図１は、インク１０２およびインク１０２と空気との間のインタフェース１０４を含む
代表的なインクジェット印刷ヘッドのノズル１００を示す。圧力パルスをインク１０２に
印加することができ、それによりインクの液滴がインタフェース１０４で形成される。圧
力パルスは圧力板、すなわちダイヤフラム、経由でインク１０２に結合される圧電（ＰＺ
Ｔ）アクチュエータに動電圧を印加することにより生成することができる。
【００２５】
印刷ヘッドのノズル１００を設計する場合計算流体力学（ＣＦＤ）コードを用いてイン
クの液滴の生成をシミュレートすることが有用である。ＣＦＤコードは流体速度、圧力、
およびインタフェース位置に関する支配的な方程式の組（すなわち二相流に関する非圧縮
性ナビエ・ストークス方程式）を解くことによりインクの液滴の生成をシミュレートする
。ナビエ・ストークス方程式は一例であり、インクの液滴の挙動を記述する他の方程式も
請求される発明の範囲を変えることなく用いることができる。
【００２６】
図２は、先行技術で実施されるような数値的シミュレーションの結果を図示する。この
例において、半径０．３の円形インク液滴が半径３、長さ２０の円筒型流路内で空気に囲
まれている。初期時間ステップｔ0において、速度は１、密度比は８７７、動粘性比は２
０８、レイノルズ数は３０．０、ウェーバー数は３１．３で、流入圧力および流出圧力は
ゼロに設定されている。このシミュレーションは均一３０×２００の四角形メッシュで実
施され、インク液滴の挙動を記述する速度、圧力、およびレベルセットを計算するために
用いられる。時間ステップはΔｔ＝０．００２となるように取られる。レベルセットは６
０時間ステップ毎に再初期設定される。図２は、初期時間ステップｔ0および後続のいく
つかの時間ステップｔ1，ｔ2，ｔ3，ｔ4，ｔ5におけるシミュレーションによる結果を示
す。このシミュレーションの結果はサイズが減少し、質量が保存されない液滴を示す。こ
のモデルに蒸発は含まれない。従ってこの液滴の縮小はシミュレーションの産物である。
本発明の一目的はこの種の産物の出現および影響を削減することにある。
【００２７】
図３は、ＣＦＤコードが発明の一実施形態においてインクの液滴が生成される物理的空
間を表すために用いることができる、境界に合わせた四辺形メッシュの例である。
【００２８】
二相流の支配方程式は連続方程式（１）、ナビエ・ストークス方程式（２）、およびレ
ベルセット対流方程式（３）を含み、これらの方程式を以下に示す。
【００２９】
【数１】
2008047109.tif 000003
【００３０】
【数２】
2008047109.tif 000004
【００３１】
【数３】
2008047109.tif 000005
これらの方程式中、ｕは速度ベクトル、ｔは時間、ρは相対密度、ｐは圧力、μは相対
動粘性、Ｒｅはレイノルズ数、Ｗｅはウェーバー数、Ｆｒはフルード数、κは曲率、δは
ディラック・デルタ関数、およびＤは変形テンソル率である。相対密度ρ、相対動粘性μ
、および曲率κは皆、レベルセット関数（φ）に関して定義される。
【００３２】
代表的な印刷ヘッドのノズル１００は軸対称である。従って、印刷ヘッドのノズル１０
０を円筒形の座標システムでモデル化すると都合が良い。結果は方位に無関係であると仮
定することが妥当である。従ってモデルは三次元から軸対称の物理的空間の二次元、Ｘ（
ｒ，ｚ）、に約分することができる。形状に合った四辺形のメッシュ３００はこの物理的
軸対称空間Ｘに作成される。
【００３３】
図４に示されるように、この物理的軸対称空間Ｘは物理的軸対称空間Ｘにおけるメッシ
ュ点を演算空間Ξ＝（ξ，η）および長方形メッシュとしてのメッシュ４００にマップす
る変換Φにより変換することができる。この演算空間Ξにおけるメッシュ４００上の支配
方程式は有限差分解析を用いて解くことができる。この変換のヤコビアン（４）および変
換マトリクス（５）は以下に示す。
【００３４】
【数４】
2008047109.tif 000006
【００３５】
【数５】
2008047109.tif 000007
【００３６】
【数６】
2008047109.tif 000008
【００３７】
【数７】
2008047109.tif 000009
円筒型座標システムにおいてｇ＝２πｒである。連続方程式（１）、ナビエ・ストーク
ス方程式（２）、およびレベルセット対流方程式（３）の数値を演算空間Ξにおいて求め
るにはヤコビアン（４）および変換マトリクス（５）を用いて方程式（１）、方程式（２
）、および方程式（３）を上記の方程式（６）に変換する必要がある。方程式（２）の第
２項は粘性項である。上記の方程式（７）は演算空間Ξで求める場合の粘性項である。
【００３８】
四辺形メッシュ上で数値アルゴリズムが公式化される。以下において、上付き記号ｎ（
またはｎ＋１）は時間ステップを意味する。量ｕn，ｐn，φnが与えられた場合、アルゴ
リズムの目的は支配方程式を満足するｕn+1、ｐn+1、およびφn+1を得ることである。発
明の一実施形態は時間的に一次の精度を有し、空間的に二次の精度を有する。発明はこの
精度を有するアルゴリズムに限定されず、他の公式を含むことができる。
【００３９】
図４に見られるように、変換Ｘ＝Φ（Ξ）は演算空間における演算メッシュとしてのメ
ッシュ４００が単位四角形、Δξ＝Δη＝１、または長方形メッシュで構成されるような
ものである。図３における境界に合わせた四辺形メッシュおよびノズル壁は単位四角形を
有する演算メッシュとしてのメッシュ４００にマップされる。
【００４０】
以下の考察を通して変数は１つ以上のメッシュのノードおよびエッジにおいて計算され
る。ノードの相対位置はカッコ内または下付き記号として（Ｚ（ｘ，ｙ）⇔Ｚx,y）変数
のペアで表される。
【００４１】
図５は、四辺形メッシュとしてのメッシュ３００の一部を図示する。図５で示す通り、
整数の時間ステップにおいて圧力ｐn（ｉ，゜゜ｊ）は四辺形メッシュとしてのメッシュ
４００のノードで計算される。さらに整数の時間ステップにおいて速度成分ｕn（ｉ，゜
ｊ）は四辺形メッシュの各セルの中心において計算される。レベルセット値φnは四辺形
メッシュとしてのメッシュ４００より奇数倍細かいメッシュ上で計算される。
【００４２】
圧力は第１メッシュ上で計算される。速度は第１メッシュよりオフセットされる第２メ
ッシュ上で計算される。第２メッシュは、第２メッシュにおける各ノードが第１メッシュ
の４ノードの中心にあるようにオフセットされる。レベルセット値は第３メッシュ上で計
算される。第３メッシュは第２メッシュおよび第１メッシュより奇数倍細かい（３、５、
７、９、等々）。第３メッシュのノードのサブセットが第２メッシュのノードに位置合わ
せされる。
【００４３】
図６（Ａ）は、均一メッシュとしてのメッシュ４００の一部を示し、ここでレベルセッ
ト値φnは３倍の微細化メッシュ上で計算される。図６（Ｂ）は、均一メッシュとしての
メッシュ４００の一部を示し、ここでレベルセット値φnは５倍の微細化メッシュ上で計
算される。図６（Ｃ）は、均一メッシュとしてのメッシュ４００の一部を示し、ここでレ
ベルセット値φnは７倍の微細化メッシュ上で計算される。
【００４４】
アルゴリズムでは変換（またはマッピング）の局部定義のみが必要であることが特記さ
れる。大域変換Ｘ＝Φ（Ξ）の存在または正確な形は重要でない。ヤコビアン（４）およ
び変換マトリクス（５）は正規のメッシュおよび奇数倍細かいメッシュ双方に対し計算さ
れる必要がある。
【００４５】
代表的なシステムで、両流体の密度および粘性は非常に差がある。インタフェース１０
４を挟んだこの急激な変化はシミュレーションにおいて問題を起こす可能性がある。従っ
てこの急激な変化を平滑化関数で置き換えることができる。平滑化はインタフェース１０
４を中心に厚み２εの空間領域で生じ得る。従ってεはインタフェース１０４のスミアリ
ングの度合いを表す。このような平滑化関数の例としては方程式（６）に示される下記の
平滑化ヘビサイド関数（８）が挙げられる。
【００４６】
【数８】
2008047109.tif 000010
代表的なεの値はメッシュ３００の四辺形セルの平均サイズに対し比例して１．７から
２．５倍であることができる。加えて、平滑化関数の導関数は方程式に示されるように平
滑化ヘビサイド関数（８）は、デルタ関数として定義することができる。
【００４７】
インタフェース１０４はレベルセット関数（φ）のゼロ・レベルセットΓで表される。
レベルセット関数（φ）はインタフェース１０４への符号付き距離として初期設定される
。すなわちレベルセット値はインク１０２側でインタフェース１０４までの最短距離で空
気側では最短距離の負となる。レベルセット関数（φ）は第３メッシュの各ノードにおい
て計算される。
【００４８】
下記の方程式（９）に示すとおりレベルセット対流方程式（φt＋ｕ・∇φ＝０）から
導かれる方程式（９）を用いてレベルセット値φn+1を更新することができる。方程式（
９）における時間中心の移流項は時間ｔnで得られるデータのみを必要とする陽公式の予
測修正スキームを用いて数値を求めることができる。
【００４９】
【数９】
2008047109.tif 000011
【００５０】
速度ベクトルを更新するには、方程式（６）から導くことのできる下記の方程式（１０
）を用いて中間値ｕ*をまず計算することができる。その後、下記の方程式（１１）を用
いて圧力フィールドｐn+1を計算することができる。下記の方程式（１２）を用い新しい
速度ベクトル・フィールドｕn+1を計算することができる。
【００５１】
【数１０】
2008047109.tif 000012
【００５２】
【数１１】
2008047109.tif 000013
【００５３】
【数１２】
2008047109.tif 000014
【００５４】
インタフェース１０４を正しく捕らえ表面張力を正確に計算するために、レベルセット
はインタフェース１０４に対する符号付き距離関数として維持する必要がある。しかしレ
ベルセットが方程式（９）により更新された場合はそうでなくなる。従ってシミュレーシ
ョンは定期的に停止され新しいレベルセット関数（φ）が再現される。
【００５５】
移流項はセル中心の予想修正スキームを用いて計算することができ、これはＪｏｈｎ  
Ｂ．ＢＥＬＬ他が「可変密度の流れに対する二次投影法」（「A Second-order Projectio
n Method for Variable-Density Flows」）、Journal of Computational Physics、101、
334-348頁、1992年、において二相流（２つの一定密度）に関し述べた非分割二次ゴドゥ
ノフ法に基づき得る。
【００５６】
下記の方程式（１３）は方程式（９）に見られる時間中心の移流項を計算するのに用い
られる修正ステップである。３倍微細化メッシュについて、定数ａおよび定数ｂは（１，
０および−１）の値を取り得る。当業者であれば方程式（１３）を修正して５、７、また
は９倍微細化メッシュなど奇数倍微細化メッシュを記述するのは自明であろう。
【００５７】
【数１３】
2008047109.tif 000015
【００５８】
図７は、均一メッシュとしてのメッシュ４００および３倍微細化メッシュの一部におけ
るレベルセット値φn、速度ｕn、時間中心のレベルセット値φn+1/2および時間中心のエ
ッジ速度ｕn+1/2の相対位置を示す。方程式（１３）は方程式（９）からの時間中心の移
流項が３倍微細化メッシュ上で計算できるかの例である。
【００５９】
2008047109.tif 000016
【００６０】
項φn+1/2（ｉ＋ａ／３＋１／６，゜ｊ＋ｂ／３）は次の方法を用いて計算される。左
側から補外すると下記の方程式（１４）が得られる。右側から補外すると下記の方程式（
１５）が得られる。
【００６１】
【数１４】
2008047109.tif 000017
【００６２】
【数１５】
2008047109.tif 000018
【００６３】
2008047109.tif 000019
【００６４】
正規導関数を得ることは中央、前方、および後方差分を下記に示した式（１６）の計算
が係わる。限定的導関数（a limiting derivative）が下記の方程式（１７）で定義され
る。二次の限定された導関数（the second-order limited ｄｅｒｉｖａｔｉｖｅ）が下
記の方程式（１８）で与えられる。接線導関数は風上差分スキームを用いて計算すること
ができる。下記の方程式（１９）は接線導関数を計算するために用いることのできる一方
法である。
【００６５】
【数１６】
2008047109.tif 000020
【００６６】
【数１７】
2008047109.tif 000021
【００６７】
【数１８】
2008047109.tif 000022
【００６８】
【数１９】
2008047109.tif 000023
【００６９】
必要な点すべてにおいて時間中心の移流流束の値を求めるのにテイラー展開を用いるの
は実際的でないかもしれない。これは図５、図６（Ａ）〜図６（Ｃ）に示されるように、
任意の時間ステップｔnにおける速度は基本メッシュのセル中心においてのみ定義される
からである。従ってテイラー展開は基本メッシュのセル・エッジにおける時間中心速度し
か与えない。代わりに下記の方程式（２０）に示すように、距離加重局部平均を用いて時
間中心の移流速度を計算する。得られた移流速度は次に下記の方程式（２１）に示すよう
に、より細かいメッシュのエッジにおける移流流束を計算するために用いられる。これら
の移流流束は下記の方程式（２２）に示すように、ゴドゥノフ型風上差分スキームで用い
てどの時間中心エッジ・レベルセット値（（１３）および（１４）より）を保持するかを
決めることができる。
【００７０】
【数２０】
2008047109.tif 000024
【００７１】
【数２１】
2008047109.tif 000025
【００７２】
【数２２】
2008047109.tif 000026
【００７３】
図８は、新しいレベルセット値φn+1がどのように奇数倍細かいメッシュ上で計算でき
るかを示すフローチャートである。図５で示すように、ある時間ステップｎにおける流体
の速度ｕnは時間ステップｎにおける流体の圧力ｐnが計算された第１メッシュからオフセ
ットされる第２メッシュ上で計算することができる。新しいレベルセット値φn+1の計算
には時間中心のエッジ速度の計算が必要かもしれない。レベルセット値が奇数倍細かいメ
ッシュ上で計算された場合、時間中心のエッジ速度もその奇数倍細かいメッシュ上で計算
される必要があるかもしれない。時間中心のエッジ速度はステップ８０２において基本メ
ッシュ上で計算することができる。これはテイラー展開およびゴドゥノフ風上差分を用い
て行なうことができる。基本メッシュ上の時間中心のエッジ速度は次に奇数倍細かいメッ
シュ上で時間中心のエッジ速度を計算するのに用いることができ、これはステップ８０４
において方程式（２０）を用い、最隣接の２から４の一次距離加重局部平均を用いて計算
することができる。次にステップ８０６において方程式（２１）を用い奇数倍細かいメッ
シュ上の時間中心のエッジ流束を計算することができる。
【００７４】
新しいレベルセット値φn+1の計算は新しいレベルセット値φn+1を計算するのに役立つ
中間値の計算が関わることができる。例えばステップ８０８において方程式（１６）を用
い、微細化メッシュ上で演算メッシュの第１軸に沿ったレベルセットの中央、前方、およ
び後方差分を計算することができる。特定ノードにおける前方差分は第１軸に沿った正の
方向における最も近い隣接ノードのレベルセット値と特定ノードにおけるレベルセット値
との間の差分である。特定ノードにおける後方差分は特定ノードにおけるレベルセット値
と第１軸に沿った負の方向における最も近い隣接ノードのレベルセット値との間の差分で
ある。中央差分は第１軸に沿った正の方向における最も近い隣接ノードと負の方向におけ
る最も近い隣接ノードのレベルセット値との差分の半分である。当業者であれば演算メッ
シュの第２軸に沿った後方、前方、および中央差分を計算するのは自明であろう。
【００７５】
前方および後方差分は、ステップ８１０において方程式（１７）を用い、微細化メッシ
ュのノードにおける演算空間の第１軸に沿ったレベルセットの限定的導関数を計算するの
に利用することができる。特定ノードにおける前方および後方差分が反対の符号を有し、
または少なくとも片方がゼロである場合、限定的導関数はゼロに等しい。それ以外の場合
、限定的導関数は前方差分または後方差分のいずれか小さい方の絶対値の２倍に等しい。
【００７６】
限定的導関数および中央差分はステップ８１２において方程式（１８）を用い、微細化
メッシュの特定ノードにおける演算空間の第１軸に沿ったレベルセットの二次の限定され
た導関数を計算するのに利用できる。第１軸は実空間の半径方向に相当することができる
。限定的導関数の小さい方または中央差分の絶対値を中央差分の符号で乗じ、二次の限定
された導関数を生み出すことができる。
【００７７】
ステップ８１２において方程式（１９）を用い、接線導関数項を計算することができる
。演算空間の第２軸は実空間の軸方向の軸に関係することができる。特定ノードにおける
第２軸に沿った速度がゼロより大きい場合、接線導関数項は第２軸に沿った速度掛けるそ
の特定ノードおよび第２軸に沿って負の方向の隣接ノードにおけるレベルセット値間の差
分に設定することができる。その特定ノードにおける第２軸に沿った速度がゼロ未満の場
合、接線導関数項は第２軸に沿った速度掛ける第２軸に沿って正の方向の隣接ノードにお
けるレベルセット値およびその特定ノードにおけるレベルセット値間の差分に設定するこ
とができる。
【００７８】
微細化メッシュの垂直エッジにおける左および右の時間中心レベルセット値はステップ
８１４において方程式（１４）および（１５）を用い計算することができる。方程式（１
４）は左からの補外で、方程式（１５）は右からの補外である。
【００７９】
微細化メッシュのエッジにおける時間中心レベルセット値はステップ８１６において方
程式（２２）を用い計算することができる。特定エッジにおける正規移流流速がゼロより
大きい場合、その特定エッジにおける時間中心レベルセット値はその特定エッジにおける
左の時間中心レベルセット値に設定できる。特定エッジにおける正規移流流速がゼロ未満
場合、その特定エッジにおける時間中心レベルセット値はその特定エッジにおける右の時
間中心レベルセット値に設定できる。特定エッジにおける速度がゼロの場合、その特定エ
ッジにおける時間中心レベルセット値はその特定エッジにおける左および右双方の時間中
心レベルセット値の平均に設定することができる。当業者であればその特定エッジに直交
するエッジにおける時間中心レベルセット値を計算するのは自明であろう。
【００８０】
演算空間における微細化メッシュ上のノードの時間中心移流項はステップ８１８におい
て方程式（１３）を用い計算することができる。方程式（１３）は時間中心移流項の高次
差分化であることができる。演算空間における微細化メッシュ上のノードの新しいレベル
セット値はステップ８２０において方程式（９）を用い計算することができる。
【００８１】
上記の値はすべて前以てまたは必要に応じて計算することができる。これらの値は微細
化メッシュのすべてのエッジに対し計算することができ、またはインタフェース近辺のエ
ッジに対してのみ計算することができる。添付クレームの精神および範囲を逸脱すること
なく、上述された方程式の高次近似を用いて上述の値を計算することができる。
【００８２】
［偶数微細化メッシュ］
図９は、レベルセット値φがセル中心の２倍微細化メッシュ上で計算される基本的演算
メッシュとしてのメッシュ４００を示す。基本メッシュにおける速度の位置は中空の円で
印され、レベルセット値φの位置は黒点で印される。いずれのレベルセット値φも速度と
同じ位置にはない。偶数微細化メッシュ上のレベルセット値φは奇数倍微細化メッシュ上
のレベルセット値φを更新するのに用いた方程式に類似した方程式を用いて更新すること
ができる。いずれの速度もレベルセット値φと同じ位置にないので、新しいレベルセット
値を計算する際使用される速度値はすべてテイラー展開を行なう前に基本メッシュ上の速
度を平均して計算する必要がある。これは数値的粘性を導入し、メッシュ解像度を削減し
、シミュレーション精度を下げる。本発明の一実施形態のように３倍微細化メッシュを用
いると、基本メッシュ上で９個中１個のレベルセット値φが速度と同じ位置になる。
【００８３】
図１０は、本発明の一実施形態で実施され得る数値シミュレーションの結果を図示する
。この例において、半径０．３の円形インク液滴が半径３、長さ２０の円筒型流路内で空
気に囲まれている。初期時間ステップｔ0において、速度は１、密度比は８７７、動粘性
比は２０８、レイノルズ数は３０．０、ウェーバー数は３１．３で、流入圧力および流出
圧力はゼロに設定されている。このシミュレーションは均一３０×２００の四角形メッシ
ュで実施され、速度と圧力を計算するのに用いられた。レベルセット向けの奇数倍微細化
メッシュは９０×６００である。時間ステップはΔｔ＝０．００２となるように取られる
。レベルセットは６０時間ステップ毎に再初期設定される。図２は、初期時間ステップｔ
0および後続のいくつかの時間ステップｔ1，ｔ2，ｔ3，ｔ4，ｔ5におけるシミュレーショ
ンによる結果を示す。このシミュレーションの結果はサイズが大幅に変わらない液滴を示
す。このように、このシミュレーションにおいては必要とする演算リソースを大幅に増や
すことなく質量の保存が大幅に増加する。
【００８４】
［システム］
発明の詳細を説明した上、本発明の１つ以上の態様を実施するために用いることのでき
る代表的なシステムを以下、図１１を参考に説明する。図１０に示されるように、システ
ムは演算リソースを提供し、コンピュータを制御する中央演算装置（以降、ＣＰＵと称す
る）１００１を含む。ＣＰＵ１００１はマイクロプロセッサ等々により実施することがで
き、グラフィックス・プロセッサおよび／あるいは数学的演算用に浮動小数点プロセッサ
を含むことができる。システム１０００はさらにシステム・メモリ１００２を含むことで
き、これはランダムアクセス・メモリ（ＲＡＭ）および読み取り専用メモリ（ＲＯＭ）の
形を取ることができる。
【００８５】
図１０に示すように、コントローラおよび周辺機器もいくつか備えることができる。入
力コントローラ１００３はキーボード、マウス、またはスタイラスなど各種入力装置１０
０４へのインタフェースを表す。スキャナ１００６と通信するスキャナ・コントローラ１
００５もあって良い。システム１０００はさらに記憶コントローラ１００７を含むことが
でき、これはオペレーティング・システム、ユーティリティ、および本発明の各種態様を
実施するプログラムの具現を含むことのできるアプリケーション向けの命令からなるプロ
グラムを記録するために用いることのできる磁気テープもしくはディスク、または光媒体
などを各々含む１つ以上の記憶装置１００８とのインタフェースとなる。さらに発明に従
い処理されたデータ、または処理されるデータを記憶するために記憶装置１００８を用い
ることができる。システム１０００はさらにブラウン管（ＣＲＴ）または薄膜トランジス
タ（ＴＦＴ）ディスプレーであり得る表示装置１０１１へのインタフェースを提供するた
めの表示コントローラ１００９を含むことができる。システム１０００はプリンタ１０１
３と通信するプリンタ・コントローラ１０１２も含むことができる。通信コントローラ１
０１４はシステム１０００がインターネット、ローカル・エリア・ネットワーク（ＬＡＮ
）、広域ネットワーク（ＷＡＮ）を含む各種ネットワークを通じ、または赤外線信号を含
む任意の適当な電磁キャリア信号を通じ遠隔装置と接続することを可能にする１つ以上の
通信装置１０１５とインタフェースすることができる。
【００８６】
図示されたシステムにおいて、主要なシステム・コンポーネントはすべて２つ以上の物
理的バスを表し得るバス１０１６に接続されることができる。しかし、各種システム・コ
ンポーネントは互いに物理的に近接していてもいなくても良い。例えば、入力データおよ
び／あるいは出力データは１つの物理的な場所から別の場所に遠隔的に伝送することがで
きる。また、発明の各種態様を実施するプログラムはネットワーク上で遠隔地（例、サー
バ）からアクセスしても良い。このようなデータおよび／あるいはプログラムは磁気テー
プもしくはディスク、または光ディスク、または送信機・受信機のペアを含む各種のコン
ピュータ読み取り可能な媒体のいずれかを通じて搬送することができる。
【００８７】
本発明はソフトウェアで都合よく実施することができる。しかし、ハードウェアによる
実施、またはソフトウェア／ハードウェアによる実施を含み他の実施ももちろん可能であ
る。ハードウェアにより実施される機能はＡＳＩＣ、デジタル信号処理回路等々を用いて
実現することができる。従ってクレームにおける「手段」の用語はソフトウェアおよびハ
ードウェアによる実施双方を対象とすることを意図している。同様に、本明細書で用いら
れる「コンピュータ読み取り可能な媒体」の用語はソフトウェア、命令のプログラムを組
み込んだハードウェア、またはその組み合わせを含む。これら実施の選択を念頭に置き、
図面およびそれらに伴う説明は当業者が必要な処理を実施するためのプログラム・コード
（すなわちソフトウェア）を書き、または回路を制作するために必要とする機能的情報を
提供するもとであると理解される。
【００８８】
発明のさらなる態様により、上述の方法またはそのステップのいずれも命令のプログラ
ム（例、ソフトウェア）に具現化することができ、命令のプログラムはコンピュータまた
は他のプロセッサ制御装置に記憶され、または搬送され、実行されることができる。ある
いは、方法またはそのステップのいずれも機能的に同等のハードウェア（例、特定用途向
集積回路（ＡＳＩＣ）、デジタル信号処理回路、等々）またはソフトウェアとハードウェ
アの組み合わせで実施することができる。
【００８９】
発明はいくつかの具体的な実施形態に関連して説明されたが、前述の説明に照らし当業
者であればさらに別の選択、変更、変形および用途が明らかになろう。このような別の選
択は例えば以下を含む。より高い次元、別の座標システム、別の流体、別の境界条件、別
の流路、別の支配的方程式、別の数値的近似、および２つより多い流体のシステムである
。このように本明細書で説明される発明は添付クレームの精神および範囲に入るこのよう
な別の選択、変更、変形および用途をすべて含むことを意図している。
【図面の簡単な説明】
【００９０】
図面において、同種の参照記号は同種の部分を指す。
【図１】発明の一実施形態がシミュレートし得る代表的なインクジェット・ヘッドのノズルを示す。
【図２】液滴の動きをシミュレートする際に質量保存がない状態を示す。
【図３】インクジェット・シミュレーション用の、境界に合わせた四辺形メッシュを示す。
【図４】四辺形の演算空間Ξにおけるメッシュ点の物理的軸対称空間Ｘへの変換Φを示す。
【図５】物理的軸対称空間Ｘにおける四辺形メッシュ上で変数が計算される相対的メッシュ点を示す。
【図６】（Ａ）、（Ｂ）、（Ｃ）は、演算空間Ξにおける四辺形メッシュ上で変数が計算される相対的メッシュ点を示す。
【図７】演算空間Ξにおける四辺形メッシュ上で変数が計算される相対的メッシュ点を示す。
【図８】液滴の動きのシミュレーションに係わるいくつかのステップを示すフローチャート。
【図９】演算空間Ξにおける四辺形メッシュ上で変数が計算される相対的メッシュ点を示す。
【図１０】本発明の一実施形態で形成された液滴の動きのシミュレーション結果を示す。
【図１１】本発明の態様を実施するために用いることができる代表的なシステムを示すブロック図である。
【符号の説明】
【００９１】
１００４…入力装置、１００６…スキャナ、１００８…記憶装置、１０１１…表示装置
、１０１３…プリンタ、１００３…入力コントローラ、１００５…スキャナ・コントロー
ラ、１００７…記憶コントローラ、１００９…表示コントローラ、１０１２…プリンタ・
コントローラ、１００１…ＣＰＵ、１０１４…通信コントローラ、１００２…システム・
メモリ、１０１５…通信装置。
【図１】
2008047109.tif 000027
【図２】
2008047109.tif 000028
【図３】
2008047109.tif 000029
【図４】
2008047109.tif 000030
【図５】
2008047109.tif 000031
【図６】
2008047109.tif 000032
【図７】
2008047109.tif 000033
【図８】
2008047109.tif 000034
【図９】
2008047109.tif 000035
【図１０】
2008047109.tif 000036
【図１１】
2008047109.tif 000037