(57)【特許請求の範囲】
【請求項１】
流体の混練状態をコンピュータで評価するための解析方法であって、
前記流体が収容されて混練される混練空間を有限個の微小な要素で分割して混練空間モデルを設定するステップ、
前記流体をモデル化した流体モデルを設定するステップ、
前記混練空間モデルに前記流体モデルを充填率１００％未満で配置するとともに必要な混練条件を設定するステップ、
前記混練条件に基づいた前記流体モデルの流動計算を行うとともに、該流体モデルに配した仮想粒子の位置情報を追跡する粒子追跡ステップ、及び、
前記仮想粒子の位置情報と、前記混練空間モデル内での前記流体モデルの理想混練状態とを比較することにより、前記流体モデルの混練度合いを計算する評価ステップを含み、
前記粒子追跡ステップの前記流動計算では、前記混練空間モデルの各要素の体積中の前記流体モデルの充填率である体積分率が定義され、
しかも、前記粒子追跡ステップは、
前記混練空間モデル内に複数個の第１仮想粒子を配置するステップ、
前記第１仮想粒子を前記流体モデルの物理量に基づいて移動させるステップ、
前記第１仮想粒子の移動先の要素の前記前記体積分率が予め定めた値以上でない場合に、この第１仮想粒子を消滅させるステップ、及び、
前記各第１仮想粒子の位置情報を記憶するステップを含み、
前記評価ステップは、前記混練空間モデル内に複数個の第２仮想粒子をランダムに固定配置するステップ、
前記粒子追跡ステップの計算結果から現在の流体モデルの存在位置を計算するステップ、
前記流体モデルの存在位置にある前記第２仮想粒子に基づいて、前記理想混練状態を計算するステップ、及び、
前記第１仮想粒子の配置と、前記第２仮想粒子の配置との一致度によって前記流体モデルの混練度合いを計算するステップを含むことを特徴とする流体の混練状態の解析方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、流体の混練状態を解析するのに役立つシミュレーション方法に関する。
【背景技術】
【０００２】
架橋前のゴムや樹脂等は、ポリマーとともに各種の添加剤や配合剤が混ぜ合わされ、バンバリーミキサー等で混練される。混練工程では、短時間で各材料が均一に混ぜ合わされることが重要である。このため、従来、ミキサーの内部で回転するローターの形状や混練空間であるチャンバー室の形状等について、種々の開発が行われているが、試作及び実験の繰り返しによって多くの開発コストが必要であった。
【０００３】
そこで、近年では、コンピュータを用いた数値シミュレーションが種々提案されている。該数値シミュレーションでは、例えば、架橋前のゴムや樹脂といった流動性を有する可塑性材料が、バンバリーミキサー内で混練されるときの流れを計算するものが提案されている（下記特許文献１参照）。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】V. Collin (1)*, E. Peuvrel-Disdier (1) ら「Numerical and Experimental Study of Dispersive Mixing of Agglomerates」
【非特許文献２】Winston Wang and Ica Manasloczower 「Analysis of Dispersive and Distributive Mixing In Terms of Minor Component Size and Spatial Distributions In Continuous Polymer Processing Equipment」
【発明の概要】
【発明が解決しようとする課題】
【０００５】
上記非特許文献２では、数値シミュレーションを行い、可塑性材料の混練状態を定量的に把握することが提案されている。すなわち、理想となる均一な混練状態と、現在の混練状態とを比較し、これらを用いて混練状態が定量化されている。
【０００６】
しかしながら、非特許文献２の方法では、理想となる均一な混練状態は、混練空間に可塑性材料が充填率１００％の状態で充填された状態で計算されている。しかしながら、現実の混練工程では、充填率は１００％未満であるから、非特許文献２で求められる理想となる混練状態は、正しい理想状態を表しているとはいえない。
【０００７】
例えば、混練空間への可塑性材料の充填率が小さくなると、充填率が１００％である完全充填状態と比較して、可塑性材料が混練空間の壁面近傍にとどまる確率が増える。換言すれば、充填率が小さくなると、可塑性材料が混練空間内で存在しやすい箇所が壁面に偏り、理想混練状態も、壁面に偏った状態になる。しかしながら、上記非特許文献２のように、理想混練状態は充填率に依存しないという仮定に基づけば、充填率が低いほど、分散状態は悪いという評価になるおそれがある。
【０００８】
本発明は、以上のような問題点に鑑み案出なされたもので、可塑性材料といった流体の混練状態を正確に評価し、例えば混練設備の開発効率等を高め得る解析方法を提供することを目的としている。
【課題を解決するための手段】
【０００９】
本発明のうち請求項１記載の発明は、流体の混練状態をコンピュータで評価するための解析方法であって、前記流体が収容されて混練される混練空間を有限個の微小な要素で分割して混練空間モデルを設定するステップ、前記流体をモデル化した流体モデルを設定するステップ、前記混練空間モデルに前記流体モデルを充填率１００％未満で配置するとともに必要な混練条件を設定するステップ、前記混練条件に基づいた前記流体モデルの流動計算を行うとともに、該流体モデルに配した仮想粒子の位置情報を追跡する粒子追跡ステップ、及び、前記仮想粒子の位置情報と、前記混練空間モデル内での前記流体モデルの理想混練状態とを比較することにより、前記流体モデルの混練度合いを計算する評価ステップを含み、前記粒子追跡ステップの前記流動計算では、前記混練空間モデルの各要素の体積中の前記流体モデルの充填率である体積分率が定義され、しかも、前記粒子追跡ステップは、前記混練空間モデル内に複数個の第１仮想粒子を配置するステップ、前記第１仮想粒子を前記流体モデルの物理量に基づいて移動させるステップ、前記第１仮想粒子の移動先の要素の前記前記体積分率が予め定めた値以上でない場合に、この第１仮想粒子を消滅させるステップ、及び、前記各第１仮想粒子の位置情報を記憶するステップを含み、前記評価ステップは、前記混練空間モデル内に複数個の第２仮想粒子をランダムに固定配置するステップ、前記粒子追跡ステップの計算結果から現在の流体モデルの存在位置を計算するステップ、前記流体モデルの存在位置にある前記第２仮想粒子に基づいて、前記理想混練状態を計算するステップ、及び、前記第１仮想粒子の配置と、前記第２仮想粒子の配置との一致度によって前記流体モデルの混練度合いを計算するステップを含むことを特徴とする。
【発明の効果】
【００１３】
本発明によれば、予め定められた混練条件に基づいた流体モデルの流動計算を行うとともに、該流体モデルに配した仮想粒子の位置情報を追跡する粒子追跡ステップ、及び、前記仮想粒子の位置情報と、前記混練空間モデル内での前記流体モデルの理想混練状態とを比較することにより、流体モデルの混練度合いを計算する評価ステップを含む。そして、評価ステップでは、前記理想混練状態を、粒子追跡ステップで計算される流体モデルの存在位置に基づいて計算することを特徴とする。このような解析方法では、流体の混練状態を正確に定量的に評価でき、混練設備等の開発を試作することなく行うことができ、開発効率を大幅に高めることができる。
【図面の簡単な説明】
【００１４】
【図１】可塑性材料を混練するバンバリーミキサーの部分断面図である。
【図２】本実施形態の処理手順の一例を示すフローチャートである。
【図３】混練空間モデルの斜視図である。
【図４】混練空間モデルの断面図である。
【図５】混練空間モデルを分解して示す断面図である。
【図６】混練空間モデル内に流体モデルと気相モデルとを混在して配置した状態を示す断面図である。
【図７】粒子追跡の処理手順の一例を示すフローチャートである。
【図８】粒子追跡を説明する粒子の線図である。
【図９】粒子追跡ステップの時間的経過を視覚化して示す線図である。
【図１０】評価ステップの一例を示すフローチャートである。
【図１１】粒子間距離を説明する線図である。
【図１２】粒子間距離のヒストグラムである。
【図１３】現実の混練状態と理想混練状態との粒子間距離のヒストグラムである。
【図１４】（ａ）は本発明の粒子追跡ステップで得られた現実の混練状態に対応する第１仮想粒子の模式図、（ｂ）はその理想混練状態に対応する第２仮想粒子の模式図、（ｃ）は従来の理想混練状態に対応する第２仮想粒子の模式図である。
【図１５】（ａ）は実施例のシミュレーション結果、（ｂ）は比較例のシミュレーション結果を示すグラフである。
【発明を実施するための形態】
【００１５】
以下、本発明の実施の一形態が図面に基づき説明される。
本発明は、コンピュータ（図示省略）を使用して、流体の混練状態を評価するための解析方法である。「混練」とは、例えば、ゴム材料や樹脂材料の成形時の前処理工程として、原材料の薬品、粉体などと液状バインダを分散させながら互いに濡らし、それらを均質にする作用乃至操作として定義される。代表的な混練工程は、図１に示されるような、バンバリーミキサー１を用いて行われる。
【００１６】
該バンバリーミキサー１は、ゴム材料等の原材料が混練される混練空間４の外周面を区画するケーシング２と、該ケーシング２内を回転する一対のロータ３、３とを具える。この混練空間４は、本実施形態では断面横向きの略８の字状の空間とされる。ただし、混練空間は、このような形状に限定して解釈されるものではない。
【００１７】
前記流体としては、安定的な流動状態とみなすことができる状態のものであれば特に限定されないが、本実施形態では、上述の通り、架橋前のゴムのように所定の粘度を有する可塑性材料が例に挙げて説明される。架橋前のゴムの場合、このような状態としては、十分に練られて約８０℃程度まで昇温した状態が相当する。なお、流体は、可塑性を有するゴムや樹脂又はエラストマー等に限定されるものではない。
【００１８】
図２には、本実施形態の処理手順の一例が示されており、以下、各ステップが順に説明される。
【００１９】
［混練空間モデル］
本実施形態のシミュレーション方法では、先ず、前記混練空間４を有限個の微小な要素（「セル」ということもある。）でモデル化した混練空間モデル５が前記コンピュータに入力、設定される（ステップＳ１）。
【００２０】
図３には、混練空間モデル５を視覚化した斜視図が、図４には、その断面図がそれぞれ示されている。混練空間モデル５は、前記ケーシング２の内周面によって規定される外周面５ｏと、前記回転する一対のロータ３、３の外周面で規定される内周面５ｉと、ロータ３の軸方向の両端側で前記外周面を閉じる両端面５ｓとで閉じられた三次元空間を有し、この空間が三次元の要素ｅによって分割されている。
【００２１】
混練空間モデル５は、オイラー要素ｅで分割（離散化）されている。要素分割は、四面体、六面体などの他、多面体セル（ポリヘドラルグリッド）といった三次元要素で行われる。そして、各要素について、流体モデルの圧力、温度及び／又は速度等の物理量が計算される。
【００２２】
混練空間モデル５は、前記外周面５ｏ及び両端面５ｓは変形しないが、内周面５ｉは、ロータ３の回転に合わせて回転することで、容積形状が変化する。本実施形態では、このような混練空間モデル５は、図５に分離させて示されるように、一対の回転可能な回転部５Ａ、５Ｂと、これらの間を継ぐ継ぎ部５Ｃと、これらが収容される外枠部５Ｄとの４つの部分に分けて構成される。
【００２３】
前記回転部５Ａ、５Ｂは、各々、円形の外周面５Ａｏ、５Ｂｏと、ロータ３の外周面に等しい前記内周面５ｉを有する筒状でモデル化されている。そして、回転部５Ａ、５Ｂは、各々、前記外枠部５Ｄに填め込まれるとともに、その中心Ｏａ、Ｏｂの周りで回転することにより、ロータ３の回転に伴う混練空間３の容積形状の変化を表現できる。
【００２４】
前記継ぎ部５Ｃは、一対の回転部５Ａ、５Ｂ間に配されて静止している。また、継ぎ部５Ｃは、各回転部５Ａ、５Ｂと接触する凹円弧面ｊを有している。継ぎ部５Ｃの凹円弧面ｊと回転部５Ａ、５Ｂの外周面５Ａｏ、５Ｂｏとは、スライディングサーフェース等の境界条件が定義され、これにより、混練空間モデル５の回転部５Ａ、５Ｂ内で生じる物理的な作用（力及び熱等）が、この凹円弧面ｊを介して継ぎ部５Ｃの中にある流体モデルとｂへと伝達される。
【００２５】
前記外枠部５Ｄは、前記回転部５Ａ、５Ｂ及び継ぎ部５Ｃを囲む筒状をなし、その軸方向両端は、前記両端面５ｓによって閉じられている。また、外枠部５Ｄも、回転部５Ａ、５Ｂ及び継ぎ部５Ｃの要素との接触面に、スライディングサーフェース等の境界条件が定義される。これにより、外枠部５Ｄ内で回転する回転部５Ａ、５Ｂ内で生じる物理的な作用（力及び熱等）が、両者の接触面を介して外枠部５Ｄにも伝えられる。なお、外枠部５Ｄは、ロータの作用によって大きなせん断が生じるため、材料の速度等をより詳細に計算すべく、回転部５Ａ、５Ｂ及び継ぎ部５Ｃよりも小さい要素で構成される。これによって、混練空間モデル５の外周面付近の流体モデルの速度プロファイルなどがより詳細に計算される。
【００２６】
［流体モデル］
次に、本実施形態では、流体モデルが設定され、コンピュータに入力される（ステップＳ２）。該流体モデルは、前記混練空間４内を流動する流体をモデル化したものであり、そのせん断粘度、比熱、熱伝導率及び比重といった物性が定義されコンピュータに入力される。流体には、解析目的に応じて種々のものが用いられるが、本実施形態では、上述の通り、可塑状態のゴムの物性が入力される。
【００２７】
前記せん断粘度は、例えば、解析対象となる可塑性材料から粘弾性特性（Ｇ'及びＧ”）が複数の温度条件で測定され、Cox-Merz則などを用いてせん断粘度に変換することで得られる。このようにして得られたせん断粘度ηは、例えば下記式（１）のべき乗法則で近似される。
η＝ｍγ'n-1  …（１）
ここで、ｍは係数、γ'はせん断速度、ｎは係数である。
【００２８】
前記比熱は、解析対象の可塑性材料から、例えば断熱型連続法（＠２５℃）にて測定され、その値が前記コンピュータに入力される。
【００２９】
前記熱伝導率は、解析対象の可塑性材料から、例えば熱線法（＠２５℃）にて測定され、その値が前記コンピュータに入力される。
【００３０】
［気相モデル］
次に、本実施形態では、気相モデルが設定され、コンピュータに入力される（ステップＳ３）。気相モデルは、混練空間内に存在する空気をモデル化したものである。すなわち、本発明では、混練空間モデル５に、流体モデルが充填率１００％未満で充填されるため、残りの空間には、空気をモデル化した気相モデルが充填される。このような気相モデルと流体モデルとを混練空間モデル５に配することにより、流体モデルの充填率が１００％未満の状態の流動計算が可能になる。
【００３１】
気相モデルには、粘度及び比重が定義される。これらは、空気の比重及び粘度に基づいて決定され、コンピュータに入力される。
【００３２】
［境界条件］
次に、本実施形態では、流動計算のシミュレーションに必要な境界条件等、各種の条件が定義される（ステップＳ４）。設定される境界条件としては、混練空間モデル５の表面での流速境界条件及び温度境界条件が挙げられる。
【００３３】
前記流速境界条件としては、下記の２通りが設定可能であり、シミュレーションの用途や精度等に応じて、いずれかが採用される。
ａ）壁面ノースリップ条件
ｂ）壁面スリップ条件
【００３４】
前記壁面ノースリップ条件では、流体モデルは、混練空間モデル５の表面での流速は常に零とされる。一方、前記壁面スリップ条件では、流体モデルは、混練空間モデル５の表面において流速を持つ。この場合、流体モデルと混練空間モデル５との接触面のスリップ現象は、例えば、慣例に従ってNavier's Lawなどを用いてシミュレートすることができる。
【００３５】
前記温度境界条件としては、次の２種類が設定可能である。
ａ）断熱条件：温度の計算を単純化するため、混練空間モデル５の表面において、熱が外に逃げない条件。
ｂ）全ての混練空間モデル５の表面温度が温調温度（例えば５０℃）に設定される条件。
【００３６】
他の条件としては、流体モデルの初期温度、混練空間モデル５の回転部５Ａ、５Ｂの回転数によって表現されるロータの回転数、混練空間モデル５の表面のスリップ率、混練空間モデル５の容積に対する流体モデルの充填率（１００％未満）などが挙げられる。さらに、他の条件としては、流動計算の初期状態、タイムステップ、内部処理でのイタレーションの反復回数、計算終了時刻などがある。また、シミュレーションにおいて、出力されるパラメータ等が決定される。
【００３７】
初期状態は、例えば、図６に示されるように、混練空間モデル５を横切る水平な境界面Ｓを基準として、それよりも上部を気相モデルの領域Ａとし、それよりも下部を流体モデルの領域Ｍとして混在配置できる。また、境界面Ｓのレベルを変えることにより、流体モデルの充填率が調節される。これらの条件は、シミュレーションの目的等に応じて任意に定められる。
【００３８】
［粒子追跡ステップ］
次に、本実施形態では、粒子追跡ステップが行われる（ステップＳ５）。粒子追跡ステップの具体的な処理の手順の一例が、図７に示されている。
【００３９】
粒子追跡ステップにおいて、コンピュータは、前記混練条件に基づいて、混練空間モデル５と流体モデルとを用いて流動計算を行う（ステップＳ５１）。該流動計算は、汎用の流体解析ソフトウエアを用いて行われる。該流動計算を行うことにより、混練空間モデル５の各要素毎に、流体の運動状態を特定する３方向（ｘ，ｙ，ｚ）の速度成分、流体の内部状態を特定する未知量である圧力ｐ、及び温度Ｔがそれぞれ計算される。これにより、流れ場が算出される。
【００４０】
次に、コンピュータは、流動計算によって算出された時刻（ｔ＝０）における流れ場を読み込む（ステップＳ５２）。この「流れ場」は、ある領域（この例では混練空間モデル）において、任意の時刻における流体の流れを特定しうる速度、圧力及び密度等の物理量が決定された場と定義される。
【００４１】
次に、コンピュータは、第１仮想粒子を混練空間モデル５内の所定の位置に所定の個数設置する（ステップＳ５３）。本実施形態において、第１仮想粒子は、数値シミュレーションにおいて、大きさ及び質量を有しない仮想の粒子として取り扱われ、流動モデルの流動計算には影響を与えないが、流動モデルの流れに従って移動するものである。従って、この第１仮想粒子の位置情報を追跡することにより、流動モデルの流動状態を調べることができる。
【００４２】
第１仮想粒子は、流体モデルの混練（分散）度合いを調べるために、好ましくは数百個以上、より好ましくは５００個以上が混練空間モデル５内に設置されるのが望ましい。第１仮想粒子の配設位置は、混練空間モデル５内であれば、任意に定められる。これらの個数や配設位置は、境界条件の設定（ステップＳ４）時に予め定められる。
【００４３】
次に、コンピュータは、各第１仮想粒子の位置における速度情報から、次の時間ステップ（ｔ＝１）での各第１仮想粒子の位置情報を計算する（ステップＳ５４）。この計算は、図８に示される通り、次のように行われる。
【００４４】
図８では一つの第１仮想粒子Ｐ１を例に挙げているが、実際には、多数設置される各第１仮想粒子について、以後の処理が行われる。まず、第１仮想粒子Ｐ１は、時刻ｔ＝０のときの位置（Ｘt，Ｙt，Ｚt）にある。また、Ｔ＝０の流れ場において、前記位置（Ｘt，Ｙt，Ｚt）での流体モデルは、各ｘ、ｙ及びｚ成分の速度情報として（Ｖx（t），Ｖy（t），Ｖz（t））を持っている。
【００４５】
コンピュータは、時間ステップ（Ｔ＝１）後の第１仮想粒子Ｐ１の位置（Ｘt+1，Ｙt+1，Ｚt+1）を下式を用いて計算する。
Ｘt+1＝Ｘt＋Ｖx（t）×Ｔ
Ｙt+1＝Ｙt＋Ｖy（t）×Ｔ
Ｚt+1＝Ｚt＋Ｖz（t）×Ｔ
【００４６】
次に、コンピュータは、前記流れ場を参照し、計算された次の時間ステップ（ｔ＝１）での第１仮想粒子Ｐ１の位置（Ｘt+1，Ｙt+1，Ｚt+1）における流体モデルの体積分率を確認する（ステップＳ５５）。
【００４７】
本発明では、混練空間モデル５内に、気相モデルと流体モデルとが混在するため、２つの流体を一度に扱う必要がある。このために、本実施形態では、自由界面の流れの計算で用いられるＶＯＦ（Volume of Fluid）法が用いられる。ＶＯＦ法では、二流体の界面の移動を直接計算するのではなく、混練空間モデル５の各要素の体積中の流体モデルの充填率である体積分率を定義して自由界面が表現される。従って、任意の要素について流体モデルの体積分率が０の場合、その要素は、全てが気相モデルで満たされてることを意味する。逆に、任意の要素について、流体モデルの体積分率が上限の１（＝１００％）の場合、その要素は、全てが流体モデルで満たされてることを意味する。
【００４８】
次に、コンピュータは、上記体積分率が予め指定した値以上か否かを判断する（ステップＳ５６）。この値には、例えば０．５程度が設定される。このステップにより、第１仮想粒子の次のタイムステップでの移動先に実質的に流体モデルが存在しているか否かを調べることができる。
【００４９】
ステップＳ５６の結果がＹＥＳの場合、コンピュータは、新たに計算された第１仮想粒子の位置における速度情報から、さらに次の時間ステップ（ｔ＝Ｎ＋１）での各第１仮想粒子の位置情報を計算する（ステップＳ５７）。
【００５０】
他方、コンピュータは、ステップＳ５６の結果がＮＯ、すなわち、第１仮想粒子の移動先の要素には、実質的に流体モデルが存在していないと判断した場合、この第１仮想粒子を消滅させる（ステップＳ５９）。つまり、気相モデルに飛散したような第１仮想粒子については、追跡を終了する。
【００５１】
次に、コンピュータは、予め定めた回数、上記処理を繰り返した否かを判断する（ステップＳ５８）。この回数は、各第１仮想粒子の追跡を開始して十分な時間が経過したか否かを基準に設定される。ステップＳ５８の結果がＮＯの場合、コンピュータは、ステップＳ５５以降を繰り返す。図９（ａ）乃至（ｆ）には、このような粒子追跡ステップを行った結果の一例として、第１仮想粒子Ｐ１の時間的経過が黒丸で示されている。図９から明らかなように、第１仮想粒子は、混練が進むにつれて分散されているのがわかる。
【００５２】
他方、ステップＳ５８の結果がＹＥＳの場合、コンピュータは、図２の評価ステップ（ステップＳ６）に戻る。評価ステップでは、コンピュータは、第１仮想粒子の位置情報と、混練空間モデル５内での流体モデルの理想混練状態とを比較することにより、流体モデルの混練度合いを計算する。この際、理想混練状態は、粒子追跡ステップ（ステップＳ５）で計算された流体モデルの存在位置に基づいて計算される。
【００５３】
図１０には、コンピュータが行う評価ステップの具体的な処理手順の一例が示されている。評価ステップでは、先ず、粒子追跡ステップで得られたデータをメモりに読み込む（ステップＳ６１）。読み込むデータとしては、各第１仮想粒子の３次元座標情報及び速度などが含まれる。
【００５４】
次に、コンピュータは、現在のステップ数Ｎ（ここでは、Ｎ＝Ｌとおく）での第１仮想粒子について、粒子間の距離を計算し、その度数分布を求める（ステップＳ６２）。
【００５５】
図１１には、任意の領域の流体モデルＭに含まれている複数の第１仮想粒子が白丸で示されている。前記ステップＳ６２では、混練空間モデル５内で追跡された全ての第１仮想粒子の組み合わせについて前記距離が計算される。また、図１２には、縦軸に度数、横軸に前記粒子間の距離を階級とした度数分布が示されている。この度数分布は、下記数１で表すことができる。
【００５６】
【数１】
0005514244.tif 000002
【００５７】
ここで、"ｐ"は、確率を表す関数、"eps"は度数分布を求めるときの階級（粒子間の距離）の階級幅の１／２、"Ｎ"は第１仮想粒子の個数、i、jは、第１仮想粒子を特定するための番号、”ｌ”は、２つの第１仮想粒子ｉ、ｊの粒子間の距離である。φijは、１又は０の変数であり、ｉ≠ｊかつ粒子が距離ｌ−epsと距離ｌ＋epsとの間にあれば１、ｉ＝ｊのときには０である。また、添え字"calc"は、現実の混練状態を示している。
【００５８】
次に、コンピュータは、上記ステップ数Ｎ（Ｎ＝Ｌ）での理想混練状態を求める（ステップＳ６３）。
【００５９】
このステップＳ６３では、まず、コンピュータは、混練空間モデル５の三次元空間の全領域を対象として、該混練空間モデル５内に複数個の第２仮想粒子をランダムに固定配置する。この第２仮想粒子も、第１仮想粒子と同様、数値シミュレーションにおいて、大きさ及び質量を有しない仮想の粒子として取り扱われ、流動モデルの流動計算には影響を与えない。この第２仮想粒子の配置は、例えば、ランダム関数などを用いて決定されても良い。また、第２仮想粒子の配設個数は、任意に定めることができるが、好ましくは第１仮想粒子と同程度配設されるのが望ましい。
【００６０】
次に、コンピュータは、前記粒子追跡ステップの計算結果から現在の流体モデルが存在している流体位置を計算する。この流体位置は、混練空間モデル５の各要素内での流体モデルの体積分率を調べ、該体積分率が予め定めた閾値（例えば５０％）以上である要素の位置として特定される。
【００６１】
次に、コンピュータは、上記ステップで得られた流体モデルが存在している流体位置にある第２仮想粒子のみに基づいて、前記理想混練状態を決定する。すなわち、第２仮想粒子のうち、前記流体位置にある第２仮想粒子を調べ、当該第２仮想粒子だけを残す。従って、流体位置以外に存在している第２仮想粒子は消滅させる。
【００６２】
以上のステップを行うことにより、コンピュータは、現在の流体モデルが存在する位置についての最適な混練状態（ランダムに分散された状態）を求めることができる。
【００６３】
次に、コンピュータは、上記ステップで計算された理想混練状態における第２仮想粒子について、粒子間の距離を計算し、その度数分布を求める（ステップＳ６４）。ステップＳ６４では、混練空間モデル５内で流体位置にある全ての第２仮想粒子の組み合わせについて前記距離が計算される。２つの第２仮想粒子ｉ、ｊの粒子間の距離は、φijで表されている。この度数分布は、下記数２で表すことができる。
【００６４】
【数２】
0005514244.tif 000003
【００６５】
ここで、数２の記号は、数１と同様であるが、第１仮想粒子を第２仮想粒子と読み替えられる。また、添え字"optimum"は、理想混練状態を意味している。
【００６６】
次に、コンピュータは、前記第１仮想粒子の位置情報と、前記流体モデルの理想混練状態、即ち、第２仮想粒子の位置情報とを比較することにより、流体モデルの混練度合いを計算する（ステップＳ６５）。
【００６７】
本実施形態において、上記ステップは、第１仮想粒子の配置と、第２仮想粒子の配置との一致度によって流体モデルの混練度合いを計算している。本実施形態では、混練度合いは、ＤＭＩで表され、下記数３によって計算することができる。
【００６８】
【数３】
0005514244.tif 000004
【００６９】
ここで、混練度合いを示すパラメータＤＭＩは、概略、図１３に示されるように、第２仮想粒子の粒子間の距離を用いて計算された理想混練状態の度数分布（確率関数ｐoptimum）と、第１仮想粒子の粒子間の距離を用いて計算された現実の混練状態の度数分布（確率関数ｐcalc）との差を求めている。従って、このパラメータＤＭＩが小さくなればなるほど、現実の混練状態が理想混練状態に近いことを意味する。
【００７０】
以上説明したように、本発明によれば、数値シミュレーションを行い、可塑性材料などの流体モデルの混練状態を求めることができる。そして、流体モデルの理想的な混練状態と現在の混練状態とを比較することにより、現在の混練状態を定量化して把握することができる。特に、本発明によれば、流体モデルの理想混練状態は、粒子追跡ステップで計算される流体モデルの存在位置に基づいて計算されるため、理想混練状態、ひいては流体の混練状態を正確に定量的に評価することができる。
【００７１】
図１４には、このような理想混練状態の従来との差異が模式化して示されている。図１４（ａ）には、ある時間ステップでの粒子追跡ステップで得られた現実の混練状態の模式図である。仮想線は、流体モデルが存在する領域を示し、白丸は、第１仮想粒子を示している。図１３（ｂ）には、図１４（ａ）と同一の時間ステップでの評価ステップで得られた理想混練状態の模式図である。仮想線は、流体モデルが存在する領域（図１４（ａ）と同じ）を示し、白丸は、第２仮想粒子を示している。本発明では、これらの第１仮想粒子及び第２仮想粒子の一致度が評価される。
【００７２】
これに対して、図１４（ｃ）は、図１４（ａ）と同一の時間ステップでの従来の理想混練状態の模式図が示されている。従来の理想混練状態では、仮想線で囲まれた流体モデルの存在位置とは無関係に第２仮想粒子が配置されており、適切ではないことがわかる。また、このような理想混練状態と現実の混練状態とを比較した結果も、信頼性の低いものとなる。
【実施例】
【００７３】
本発明の効果を確認するために、図３乃至５の混練空間モデルを用いて可塑性のゴム材料を流体モデルとした混練シミュレーションが行われた。主なシミュレーション条件は、次の通りである。
流体モデルの充填率：７０％及び８０％の２種類（体積分率）
ローターに対応している回転部の回転数：３０rpm
第１、第２仮想粒子の各個数１００００個
第１、第２仮想粒子の初期配設位置：混練空間モデルの中央位置
混練時間：実時間で２０分間
【００７４】
テストの結果は、図１５に示されている。図１５（ａ）は本発明方法を適用した実施例、図１５（ｂ）は理想混練状態を流体モデルの存在位置と無関係に求めた従来方法を適用した比較例、それぞれについてのＤＭＩの時間の変化を示している。比較例では、充填率が７０％と低い場合、ＤＭＩが悪い（値が大きい）傾向が強く現れている。一方、実施例では、充填率によるＤＭＩの差が小さいこと、かつ、充填率が７０％と低い場合でも、ＤＭＩが悪く見積もられるという不具合が改善されていることが確認できた。
【符号の説明】
【００７５】
１  バンバリーミキサー
２  ケーシング
３  ロータ
４  混練空間
５  混練空間モデル
【図１】
0005514244.tif 000005
【図２】
0005514244.tif 000006
【図５】
0005514244.tif 000007
【図６】
0005514244.tif 000008
【図７】
0005514244.tif 000009
【図８】
0005514244.tif 000010
【図９】
0005514244.tif 000011
【図１０】
0005514244.tif 000012
【図１１】
0005514244.tif 000013
【図１２】
0005514244.tif 000014
【図１３】
0005514244.tif 000015
【図１４】
0005514244.tif 000016
【図１５】
0005514244.tif 000017
【図３】
0005514244.tif 000018
【図４】
0005514244.tif 000019