計算装置、挙動シミュレーション方法及びそのプログラム
【課題】容器内の粒子の挙動シミュレーションを精度よく行うことができる計算装置、挙動シミュレーション方法及びそのプログラムを提供する。
【解決手段】挙動シミュレーション方法に用いられる数値計算対象領域１は、流体４と、流体４中に分散する粒子５とを備え、粒子５の大きさより小さい空間要素２で分割されている。挙動シミュレーション方法は、現在時刻Ｔに算出された粒子５が受ける力を粒子５に対応する空間要素２に設定し、Ｔ〜Ｔ＋Δｔの間の流体４’（粒子５と流体４を一体とみなした流体）の流体計算を行って、Ｔ＋Δｔにおける粒子５の速度及び位置を算出・出力する。その後、現在時刻Ｔの値をＴ＋Δｔの値に更新して、上述した数値計算処理を繰り返す。
【解決手段】挙動シミュレーション方法に用いられる数値計算対象領域１は、流体４と、流体４中に分散する粒子５とを備え、粒子５の大きさより小さい空間要素２で分割されている。挙動シミュレーション方法は、現在時刻Ｔに算出された粒子５が受ける力を粒子５に対応する空間要素２に設定し、Ｔ〜Ｔ＋Δｔの間の流体４’（粒子５と流体４を一体とみなした流体）の流体計算を行って、Ｔ＋Δｔにおける粒子５の速度及び位置を算出・出力する。その後、現在時刻Ｔの値をＴ＋Δｔの値に更新して、上述した数値計算処理を繰り返す。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、計算装置、挙動シミュレーション方法及びそのプログラムに関し、特に、容器内の粒子の運動状態を数値計算により求める計算装置、挙動シミュレーション方法及びそのプログラムに関する。
【背景技術】
【０００２】
従来より、印刷トナー容器内粒子の挙動シミュレーションを行う計算装置が知られている。具体的には、粒子粒子間、粒子攪拌機関の応力から粒子の混合状態について数値計算を行い、容器内粒子が受ける力を帯電量と粒子位置から求め、時間発展を求める。すなわち、いわゆる差分法により印刷トナー容器内粒子の挙動シミュレーションを行う方法が開示されている（例えば、特許文献１参照）。
【特許文献１】特開２０００−２１４１３４号公報
【発明の開示】
【発明が解決しようとする課題】
【０００３】
しかしながら、従来の計算装置は、印刷トナーが乾式トナーであるときのシミュレーションを行うものであるため、容器内粒子の運動に流体による効果が入らないという問題がある。特に、容器内粒子が受ける力の計算に粒子速度が用いられていないため、容器内粒子の運動に流体の効果が無視できない系においては、容器内粒子の挙動シミュレーションの精度が悪くなる傾向がある。
【０００４】
本発明の目的は、容器内の粒子の挙動シミュレーションを精度よく行うことができる計算装置、挙動シミュレーション方法及びそのプログラムを提供することにある。
【課題を解決するための手段】
【０００５】
上記目的を達成するために、請求項１記載の計算装置は、容器内の流体中に分散する粒子の挙動シミュレーションを行う計算装置であって、前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割手段と、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出手段と、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定手段と、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算手段と、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出手段と、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出手段と、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出手段とを備えることを特徴とする。
【０００６】
上記目的を達成するために、請求項７記載の挙動シミュレーション方法は、容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法であって、前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割ステップと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出ステップと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定ステップと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算ステップと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出ステップと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出ステップと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出ステップとを備えることを特徴とする。
【０００７】
上記目的を達成するために、請求項８記載のプログラムは、容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法をコンピュータにより実行させるプログラムであって、前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割モジュールと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出モジュールと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定モジュールと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算モジュールと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出モジュールと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出モジュールと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出モジュールとを備えることを特徴とする。
【発明の効果】
【０００８】
本発明によれば、容器内の流体中に分散する粒子の挙動シミュレーションを行う際に、容器内の領域を数値計算対象領域とし、この数値計算対象領域を粒子よりも小さい空間要素に分割し、現在時刻から微小時間経過までの間の粒子の運動状態の数値計算を行う際に、現在時刻における粒子の位置・速度に基づいて、その粒子が受ける力を算出し、その算出された粒子が受ける力を、その粒子に対応する上記空間要素に設定し、粒子及び流体を一体の流体として現在時刻から微小時間の間の流体計算を行い、微小時間後の上記一体の流体の速度分布を算出し、その算出された速度分布のうち、粒子に対応する空間要素の速度を加重平均して現在時刻から微小時間経過後の粒子の速度を算出し、算出された粒子速度から上記微小時間経過後の粒子の位置を算出するので、容器内の粒子の挙動シミュレーションを精度よく行うことができる。
【発明を実施するための最良の形態】
【０００９】
以下、本発明の実施の形態を図面を用いて詳述する。
【００１０】
図１は、本実施の形態に係る挙動シミュレーション方法に用いられる数値計算対象領域を概略的に示す図である。
【００１１】
図１において、数値計算対象領域１は、容器内の流体４中に分散する粒子５の挙動シミュレーションを行う、図４につき後述する計算装置としてのコンピュータ１００において、差分法を用いて粒子５の運動状態の時間発展を行うのに用いられる領域である。具体的には、数値計算対象領域１は、流体４と、流体４中に分散する粒子５と、流体４及び粒子５の存在しない領域３とを備える。この数値計算対象領域１の大きさは流体４及び粒子５を含む容器の大きさを示すものであり、ユーザの入力データに基づいて設定される。また、流体４及び粒子５の物理量の値もユーザの入力データに基づいて設定される。
【００１２】
また、数値計算対象領域１は、粒子５の大きさより小さい空間要素２で分割されている。この空間要素２の大きさも、ユーザの入力データに基づいて設定される。
【００１３】
後述する図３の処理で算出される粒子５が受ける力は、粒子５に対応する空間要素２に設定する。ここで、空間要素２の１つである領域７のように、粒子５が空間要素２の全てを覆っていれば粒子５が受ける力をその空間要素２に設定すればよい。しかし、空間要素２の１つである領域６のように、粒子５が空間要素２の一部を覆っている場合はその割合に応じて粒子５が受ける力を設定する。これにより、粒子５を流体と見做すことができ、図２に示すように、数値計算対象領域１中の粒子５が移動をしているときに、粒子５と流体４と一体として流体計算を行うことができる。
【００１４】
また、各空間要素２に設定する力はこの流体計算で用いられる式に応じて必要があれば粒子体積で割った値を設定する。さらに、上記流体計算の式に応じて必要があれば粒子５に対応する空間要素２に粒子５に対応する物理量を設定する。この場合設定する物理量としては、粒子５の質量に応じた値が挙げられる。
【００１５】
このとき同時に、流体４に対応する空間要素２に、流体４に対応する物理量（例えば流体４の質量）を設定してもよい。この場合、粒子５に対応する領域の粘性・質量密度を流体４のそれらとは異なる値に設定するのが望ましい。
【００１６】
尚、図２では、粒子５のみが移動し、流体４自体は移動しない場合を示しているが、容器内で粒子５のみならず、流体４も移動する点を粒子５の挙動シミュレーションを行うにあたり考慮するようにしてもよい。
【００１７】
図３は、図１における粒子５の運動状態計算処理の手順を示すフローチャートである。
【００１８】
図３において、まず、微小時間Δｔ、粒子５の物理量及び流体４の物理量、総計算時間Ｔａｌｌ等の入力データを読み込み（ステップＳ１００）、現在時刻Ｔ等の変数の初期化を行う（ステップＳ１０１）。本実施の形態では、現在時刻Ｔの初期値は０とする。
【００１９】
次に、現在時刻Ｔでの粒子５が受ける力の計算を、現在時刻Ｔにおける粒子５の位置及び速度に基づいて行う（ステップＳ１０２）。ここで、現在時刻Ｔにおける粒子５の位置及び速度としては、本処理開始直後の処理ではステップＳ１００でユーザ入力された値又はステップＳ１０１で初期化された値が用いられる。しかし、それ以降の処理では、このフローチャートの前回のループにおいて後述するステップＳ１０８で出力された粒子５の位置及び速度の値が用いられる。
【００２０】
その後、粒子５が受ける力を粒子５に対応する空間要素２に設定する（ステップＳ１０３）。これにより、粒子５を流体と見做すことができ、この粒子５と流体４とを一体の流体（以下「流体４’」という）と仮定して、以下のステップＳ１０４における数値計算対象領域１の流体計算を行うことができる。
【００２１】
ステップＳ１０３の設定が終了した後、Ｔ〜Ｔ＋Δｔの間の流体計算を数値計算対象領域１の流体４’について行う（ステップＳ１０４）。ここで用いられる流体計算方法はＭＡＣ法・ＳＩＭＰＬＥ法など流体の数値解析手法として公知の方法を用いればよい。
【００２２】
その後、Ｔ＋Δｔにおける流体４’の速度分布を算出する（ステップＳ１０５）。また、算出された速度分布のうち、粒子５に対応する空間要素２の速度分布を加重平均し、Ｔ＋Δｔにおける粒子５の速度を算出する（ステップＳ１０６）。さらに、その粒子速度からＴ＋Δｔにおける粒子５の位置を算出する（ステップＳ１０７）。
【００２３】
上記算出された粒子５の速度及び位置のデータを出力する出力処理を行い（ステップＳ１０８）、現在時刻Ｔの値をＴ＋Δｔの値に更新する（ステップＳ１０９）。
【００２４】
上記更新後の現在時刻Ｔが、ステップＳ１００で読み込んだ総計算時間Ｔａｌｌの値未満であるときは（ステップＳ１１０でＮＯ）、ステップＳ１０２からの処理を繰り返す。一方、上記総計算時間Ｔａｌｌの値以上であるときは（ステップＳ１１０でＹＥＳ）、そのまま本処理を終了する。
【００２５】
本処理によれば、現在時刻Ｔに算出された粒子５が受ける力を粒子５に対応する空間要素２に設定し（ステップＳ１０３）、Ｔ〜Ｔ＋Δｔの間の流体４’（粒子５と流体４を一体とみなした流体）の流体計算を行う（ステップＳ１０４）。その後、Ｔ＋Δｔにおける粒子５の速度及び位置を算出・出力する（ステップＳ１０５〜Ｓ１０８）。その後、現在時刻Ｔの値をＴ＋Δｔの値に更新して(ステップＳ１０９)、ステップＳ１０３〜Ｓ１０８の数値計算処理を繰り返す。これにより、容器内の粒子、より具体的には、数値計算対象領域１にある流体４中の粒子５の挙動シミュレーションを精度よく行うことができる。
【００２６】
図３のステップＳ１０２で粒子５が受ける力として計算される値としては、数値計算対象領域１中の流体４の加速度と粒子５の質量の積に依存した力や、粒子５の電荷及び電場の積に依存した力が挙げられる。さらには、粒子５と別の粒子との距離と相対速度に依存した力や、及び粒子５と計算領域中に設定された固体領域との距離と粒子５の速度に依存した力が挙げられる。本実施の形態では、これら４つのうちの少なくとも１つが含まれるように計算される。
【００２７】
ここで、数値計算対象領域１の加速度とは、数値計算対象領域１中で流体４’が移動するように設定した場合に、その移動に伴う加速度である。例えば、重力の影響を考慮して上記計算処理を実行する場合は、その加速度は重力加速度となる。また、粒子５の帯電を考慮して上記計算処理を実行する場合は、数値計算対象領域１には電場が形成されるので、帯電した粒子５が受けるクーロン力（粒子５の電荷及び数値計算対象領域１の電場の積）を粒子５が受ける力に含むように設定するのが望ましい。
【００２８】
また、容器内の粒子間の影響を考慮して上記計算処理を実行する場合は、一方の粒子が受ける力に、他方の粒子の位置と速度に依存する力を含むように設定するのが好ましい。例えば２つの粒子が互いに接触している場合、それぞれの粒子が離れる方向の力を含むのが望ましい。また、例えば、２つの粒子が、それらの位置と速度から微小時間中に接触することが予測される場合、一方の粒子が受ける力として、他方の粒子を離す方向の力を含むのが望ましい。
【００２９】
このように設定することにより、流体中の粒子の挙動シミュレーションをより高精度に行うことができる。
【００３０】
同様に数値計算対象領域内に固体領域を設定する場合、その固体領域と粒子の距離及び速度から粒子が受ける力を求めるのが望ましい。これにより、例えば容器内に粒子攪拌機関等が設置されている場合の流体中の粒子の挙動シミュレーションをより高精度に行うことができる。
【００３１】
また、粒子が受ける力を式１により求めても良い。ここで、Ｇｉはｉ番目の粒子が受ける力である。Ｆｉはｉ番目の粒子が受ける力の関数であり粒子の位置ｒの関数である。ｒｊは式２により求める。関数Ｆは粒子の位置ｒの関数であれば特に限定しない。Ｆｉ（ｒｊ）の値は正でもよいし負でも良い。これにより粒子相互の引力・斥力相互作用を考慮した数値計算が可能となる。
【００３２】
【数１】
【００３３】
【数２】
【００３４】
ここで、ｒ０はｉ番目の粒子の位置、νｉはｉ番目の粒子の速度、ｄｔは微小時間のｎ分の１の値、ｍｉはｉ番目の粒子の質量、νｊはｊ番目の粒子の速度、ν’は粒子速度の時間微分値である。
【００３５】
図４は、本発明の実施の形態に係る計算装置としてのコンピュータの構成を概略的に示すブロック図である。
【００３６】
本実施の形態に係る計算装置は、図４に示すコンピュータ１００により実現される。図４において、ＣＰＵ１０１は、プログラムに基づいて各種処理を実行し、装置各部を制御する。表示装置１０２は処理結果を表示する。入力装置１０３は、キーボードやマウスなどを備え、コマンドやデータを入力する。１次記憶装置１０４は、ＲＡＭなどで構成され、ワークメモリとして利用される。２次記憶装置１０５は、ＨＤやＦＤなどの不揮発記憶媒体を用いてデータやプログラムを記憶する。通信装置１０６は、専用通信回線やインターネットなどを介して他の装置とデータを送受信する。内部バス１０７は上述した装置各部を接続する。
【００３７】
２次記憶装置１０５には、図３のフローチャートにつき説明した処理手順を実現するプログラムが格納されている。このプログラムは、通信装置１０６から受信するようにしてもよい。このプログラムは、２次記憶装置１０５又は通信装置１０６から１次記憶装置１０４にロードされ、ＣＰＵ１０１により実行される。また、処理対象となるデータも、２次記憶装置１０５又は通信装置１０６から読み込まれるものとする。
【００３８】
なお、本発明の目的は、上記実施形態の機能を実現するソフトウェアのプログラムコードを記録した記憶媒体を、システム或いは装置に供給し、そのシステム或いは装置のコンピュータ（またはＣＰＵやＭＰＵ等）が記憶媒体に格納されたプログラムコードを読み出して実行することによっても達成される。
【００３９】
この場合、記憶媒体から読み出されたプログラムコード自体が前述した実施の形態の機能を実現することになり、そのプログラムコードおよび該プログラムコードを記憶した記憶媒体は本発明を構成することになる。
【００４０】
プログラムコードを供給するための記憶媒体としては、例えば、ＲＯＭ、フロッピー（登録商標）ディスク、ＰＣＭＣＩＡカードやコンパクトフラッシュ（登録商標）等のメモリカード、ハードディスク、マイクロＤＡＴ、光磁気ディスク、ＣＤ−ＲやＣＤ−ＲＷ等の光ディスク、ＤＶＤ等の相変化型光ディスク等で構成されてもよい。また、プログラムコードを、ネットワークを介してダウンロードしてもよい。
【００４１】
また、コンピュータが読み出したプログラムコードを実行することにより、上記実施の形態の機能が実現されるだけでなく、そのプログラムコードの指示に基づき、コンピュータ上で稼動しているＯＳ（オペレーティングシステム）等が実際の処理の一部または全部を行い、その処理によって前述した実施形態の機能が実現される場合も含まれることはいうまでもない。
【００４２】
さらに、記憶媒体から読み出されたプログラムコードが、コンピュータに挿入された機能拡張ボードやコンピュータに接続された機能拡張ユニットに備わるメモリに書き込まれた後、そのプログラムコードの指示に基づき、その機能拡張ボードや機能拡張ユニットに備わるＣＰＵ等が実際の処理の一部または全部を行い、その処理によって前述した実施形態の機能が実現される場合も含まれる。
【図面の簡単な説明】
【００４３】
【図１】本実施の形態に係る挙動シミュレーション方法に用いられる数値計算対象領域を概略的に示す図である。
【図２】図１の数値計算対象領域中の粒子が移動する状態を概略的に示す図である。
【図３】図１における粒子の運動状態計算処理の手順を示すフローチャートである。
【図４】本発明の実施の形態に係る計算装置としてのコンピュータの構成を概略的に示すブロック図である。
【符号の説明】
【００４４】
１ 数値計算対象領域
２ 空間要素
４ 流体
５ 粒子
【技術分野】
【０００１】
本発明は、計算装置、挙動シミュレーション方法及びそのプログラムに関し、特に、容器内の粒子の運動状態を数値計算により求める計算装置、挙動シミュレーション方法及びそのプログラムに関する。
【背景技術】
【０００２】
従来より、印刷トナー容器内粒子の挙動シミュレーションを行う計算装置が知られている。具体的には、粒子粒子間、粒子攪拌機関の応力から粒子の混合状態について数値計算を行い、容器内粒子が受ける力を帯電量と粒子位置から求め、時間発展を求める。すなわち、いわゆる差分法により印刷トナー容器内粒子の挙動シミュレーションを行う方法が開示されている（例えば、特許文献１参照）。
【特許文献１】特開２０００−２１４１３４号公報
【発明の開示】
【発明が解決しようとする課題】
【０００３】
しかしながら、従来の計算装置は、印刷トナーが乾式トナーであるときのシミュレーションを行うものであるため、容器内粒子の運動に流体による効果が入らないという問題がある。特に、容器内粒子が受ける力の計算に粒子速度が用いられていないため、容器内粒子の運動に流体の効果が無視できない系においては、容器内粒子の挙動シミュレーションの精度が悪くなる傾向がある。
【０００４】
本発明の目的は、容器内の粒子の挙動シミュレーションを精度よく行うことができる計算装置、挙動シミュレーション方法及びそのプログラムを提供することにある。
【課題を解決するための手段】
【０００５】
上記目的を達成するために、請求項１記載の計算装置は、容器内の流体中に分散する粒子の挙動シミュレーションを行う計算装置であって、前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割手段と、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出手段と、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定手段と、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算手段と、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出手段と、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出手段と、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出手段とを備えることを特徴とする。
【０００６】
上記目的を達成するために、請求項７記載の挙動シミュレーション方法は、容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法であって、前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割ステップと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出ステップと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定ステップと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算ステップと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出ステップと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出ステップと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出ステップとを備えることを特徴とする。
【０００７】
上記目的を達成するために、請求項８記載のプログラムは、容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法をコンピュータにより実行させるプログラムであって、前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割モジュールと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出モジュールと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定モジュールと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算モジュールと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出モジュールと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出モジュールと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出モジュールとを備えることを特徴とする。
【発明の効果】
【０００８】
本発明によれば、容器内の流体中に分散する粒子の挙動シミュレーションを行う際に、容器内の領域を数値計算対象領域とし、この数値計算対象領域を粒子よりも小さい空間要素に分割し、現在時刻から微小時間経過までの間の粒子の運動状態の数値計算を行う際に、現在時刻における粒子の位置・速度に基づいて、その粒子が受ける力を算出し、その算出された粒子が受ける力を、その粒子に対応する上記空間要素に設定し、粒子及び流体を一体の流体として現在時刻から微小時間の間の流体計算を行い、微小時間後の上記一体の流体の速度分布を算出し、その算出された速度分布のうち、粒子に対応する空間要素の速度を加重平均して現在時刻から微小時間経過後の粒子の速度を算出し、算出された粒子速度から上記微小時間経過後の粒子の位置を算出するので、容器内の粒子の挙動シミュレーションを精度よく行うことができる。
【発明を実施するための最良の形態】
【０００９】
以下、本発明の実施の形態を図面を用いて詳述する。
【００１０】
図１は、本実施の形態に係る挙動シミュレーション方法に用いられる数値計算対象領域を概略的に示す図である。
【００１１】
図１において、数値計算対象領域１は、容器内の流体４中に分散する粒子５の挙動シミュレーションを行う、図４につき後述する計算装置としてのコンピュータ１００において、差分法を用いて粒子５の運動状態の時間発展を行うのに用いられる領域である。具体的には、数値計算対象領域１は、流体４と、流体４中に分散する粒子５と、流体４及び粒子５の存在しない領域３とを備える。この数値計算対象領域１の大きさは流体４及び粒子５を含む容器の大きさを示すものであり、ユーザの入力データに基づいて設定される。また、流体４及び粒子５の物理量の値もユーザの入力データに基づいて設定される。
【００１２】
また、数値計算対象領域１は、粒子５の大きさより小さい空間要素２で分割されている。この空間要素２の大きさも、ユーザの入力データに基づいて設定される。
【００１３】
後述する図３の処理で算出される粒子５が受ける力は、粒子５に対応する空間要素２に設定する。ここで、空間要素２の１つである領域７のように、粒子５が空間要素２の全てを覆っていれば粒子５が受ける力をその空間要素２に設定すればよい。しかし、空間要素２の１つである領域６のように、粒子５が空間要素２の一部を覆っている場合はその割合に応じて粒子５が受ける力を設定する。これにより、粒子５を流体と見做すことができ、図２に示すように、数値計算対象領域１中の粒子５が移動をしているときに、粒子５と流体４と一体として流体計算を行うことができる。
【００１４】
また、各空間要素２に設定する力はこの流体計算で用いられる式に応じて必要があれば粒子体積で割った値を設定する。さらに、上記流体計算の式に応じて必要があれば粒子５に対応する空間要素２に粒子５に対応する物理量を設定する。この場合設定する物理量としては、粒子５の質量に応じた値が挙げられる。
【００１５】
このとき同時に、流体４に対応する空間要素２に、流体４に対応する物理量（例えば流体４の質量）を設定してもよい。この場合、粒子５に対応する領域の粘性・質量密度を流体４のそれらとは異なる値に設定するのが望ましい。
【００１６】
尚、図２では、粒子５のみが移動し、流体４自体は移動しない場合を示しているが、容器内で粒子５のみならず、流体４も移動する点を粒子５の挙動シミュレーションを行うにあたり考慮するようにしてもよい。
【００１７】
図３は、図１における粒子５の運動状態計算処理の手順を示すフローチャートである。
【００１８】
図３において、まず、微小時間Δｔ、粒子５の物理量及び流体４の物理量、総計算時間Ｔａｌｌ等の入力データを読み込み（ステップＳ１００）、現在時刻Ｔ等の変数の初期化を行う（ステップＳ１０１）。本実施の形態では、現在時刻Ｔの初期値は０とする。
【００１９】
次に、現在時刻Ｔでの粒子５が受ける力の計算を、現在時刻Ｔにおける粒子５の位置及び速度に基づいて行う（ステップＳ１０２）。ここで、現在時刻Ｔにおける粒子５の位置及び速度としては、本処理開始直後の処理ではステップＳ１００でユーザ入力された値又はステップＳ１０１で初期化された値が用いられる。しかし、それ以降の処理では、このフローチャートの前回のループにおいて後述するステップＳ１０８で出力された粒子５の位置及び速度の値が用いられる。
【００２０】
その後、粒子５が受ける力を粒子５に対応する空間要素２に設定する（ステップＳ１０３）。これにより、粒子５を流体と見做すことができ、この粒子５と流体４とを一体の流体（以下「流体４’」という）と仮定して、以下のステップＳ１０４における数値計算対象領域１の流体計算を行うことができる。
【００２１】
ステップＳ１０３の設定が終了した後、Ｔ〜Ｔ＋Δｔの間の流体計算を数値計算対象領域１の流体４’について行う（ステップＳ１０４）。ここで用いられる流体計算方法はＭＡＣ法・ＳＩＭＰＬＥ法など流体の数値解析手法として公知の方法を用いればよい。
【００２２】
その後、Ｔ＋Δｔにおける流体４’の速度分布を算出する（ステップＳ１０５）。また、算出された速度分布のうち、粒子５に対応する空間要素２の速度分布を加重平均し、Ｔ＋Δｔにおける粒子５の速度を算出する（ステップＳ１０６）。さらに、その粒子速度からＴ＋Δｔにおける粒子５の位置を算出する（ステップＳ１０７）。
【００２３】
上記算出された粒子５の速度及び位置のデータを出力する出力処理を行い（ステップＳ１０８）、現在時刻Ｔの値をＴ＋Δｔの値に更新する（ステップＳ１０９）。
【００２４】
上記更新後の現在時刻Ｔが、ステップＳ１００で読み込んだ総計算時間Ｔａｌｌの値未満であるときは（ステップＳ１１０でＮＯ）、ステップＳ１０２からの処理を繰り返す。一方、上記総計算時間Ｔａｌｌの値以上であるときは（ステップＳ１１０でＹＥＳ）、そのまま本処理を終了する。
【００２５】
本処理によれば、現在時刻Ｔに算出された粒子５が受ける力を粒子５に対応する空間要素２に設定し（ステップＳ１０３）、Ｔ〜Ｔ＋Δｔの間の流体４’（粒子５と流体４を一体とみなした流体）の流体計算を行う（ステップＳ１０４）。その後、Ｔ＋Δｔにおける粒子５の速度及び位置を算出・出力する（ステップＳ１０５〜Ｓ１０８）。その後、現在時刻Ｔの値をＴ＋Δｔの値に更新して(ステップＳ１０９)、ステップＳ１０３〜Ｓ１０８の数値計算処理を繰り返す。これにより、容器内の粒子、より具体的には、数値計算対象領域１にある流体４中の粒子５の挙動シミュレーションを精度よく行うことができる。
【００２６】
図３のステップＳ１０２で粒子５が受ける力として計算される値としては、数値計算対象領域１中の流体４の加速度と粒子５の質量の積に依存した力や、粒子５の電荷及び電場の積に依存した力が挙げられる。さらには、粒子５と別の粒子との距離と相対速度に依存した力や、及び粒子５と計算領域中に設定された固体領域との距離と粒子５の速度に依存した力が挙げられる。本実施の形態では、これら４つのうちの少なくとも１つが含まれるように計算される。
【００２７】
ここで、数値計算対象領域１の加速度とは、数値計算対象領域１中で流体４’が移動するように設定した場合に、その移動に伴う加速度である。例えば、重力の影響を考慮して上記計算処理を実行する場合は、その加速度は重力加速度となる。また、粒子５の帯電を考慮して上記計算処理を実行する場合は、数値計算対象領域１には電場が形成されるので、帯電した粒子５が受けるクーロン力（粒子５の電荷及び数値計算対象領域１の電場の積）を粒子５が受ける力に含むように設定するのが望ましい。
【００２８】
また、容器内の粒子間の影響を考慮して上記計算処理を実行する場合は、一方の粒子が受ける力に、他方の粒子の位置と速度に依存する力を含むように設定するのが好ましい。例えば２つの粒子が互いに接触している場合、それぞれの粒子が離れる方向の力を含むのが望ましい。また、例えば、２つの粒子が、それらの位置と速度から微小時間中に接触することが予測される場合、一方の粒子が受ける力として、他方の粒子を離す方向の力を含むのが望ましい。
【００２９】
このように設定することにより、流体中の粒子の挙動シミュレーションをより高精度に行うことができる。
【００３０】
同様に数値計算対象領域内に固体領域を設定する場合、その固体領域と粒子の距離及び速度から粒子が受ける力を求めるのが望ましい。これにより、例えば容器内に粒子攪拌機関等が設置されている場合の流体中の粒子の挙動シミュレーションをより高精度に行うことができる。
【００３１】
また、粒子が受ける力を式１により求めても良い。ここで、Ｇｉはｉ番目の粒子が受ける力である。Ｆｉはｉ番目の粒子が受ける力の関数であり粒子の位置ｒの関数である。ｒｊは式２により求める。関数Ｆは粒子の位置ｒの関数であれば特に限定しない。Ｆｉ（ｒｊ）の値は正でもよいし負でも良い。これにより粒子相互の引力・斥力相互作用を考慮した数値計算が可能となる。
【００３２】
【数１】
【００３３】
【数２】
【００３４】
ここで、ｒ０はｉ番目の粒子の位置、νｉはｉ番目の粒子の速度、ｄｔは微小時間のｎ分の１の値、ｍｉはｉ番目の粒子の質量、νｊはｊ番目の粒子の速度、ν’は粒子速度の時間微分値である。
【００３５】
図４は、本発明の実施の形態に係る計算装置としてのコンピュータの構成を概略的に示すブロック図である。
【００３６】
本実施の形態に係る計算装置は、図４に示すコンピュータ１００により実現される。図４において、ＣＰＵ１０１は、プログラムに基づいて各種処理を実行し、装置各部を制御する。表示装置１０２は処理結果を表示する。入力装置１０３は、キーボードやマウスなどを備え、コマンドやデータを入力する。１次記憶装置１０４は、ＲＡＭなどで構成され、ワークメモリとして利用される。２次記憶装置１０５は、ＨＤやＦＤなどの不揮発記憶媒体を用いてデータやプログラムを記憶する。通信装置１０６は、専用通信回線やインターネットなどを介して他の装置とデータを送受信する。内部バス１０７は上述した装置各部を接続する。
【００３７】
２次記憶装置１０５には、図３のフローチャートにつき説明した処理手順を実現するプログラムが格納されている。このプログラムは、通信装置１０６から受信するようにしてもよい。このプログラムは、２次記憶装置１０５又は通信装置１０６から１次記憶装置１０４にロードされ、ＣＰＵ１０１により実行される。また、処理対象となるデータも、２次記憶装置１０５又は通信装置１０６から読み込まれるものとする。
【００３８】
なお、本発明の目的は、上記実施形態の機能を実現するソフトウェアのプログラムコードを記録した記憶媒体を、システム或いは装置に供給し、そのシステム或いは装置のコンピュータ（またはＣＰＵやＭＰＵ等）が記憶媒体に格納されたプログラムコードを読み出して実行することによっても達成される。
【００３９】
この場合、記憶媒体から読み出されたプログラムコード自体が前述した実施の形態の機能を実現することになり、そのプログラムコードおよび該プログラムコードを記憶した記憶媒体は本発明を構成することになる。
【００４０】
プログラムコードを供給するための記憶媒体としては、例えば、ＲＯＭ、フロッピー（登録商標）ディスク、ＰＣＭＣＩＡカードやコンパクトフラッシュ（登録商標）等のメモリカード、ハードディスク、マイクロＤＡＴ、光磁気ディスク、ＣＤ−ＲやＣＤ−ＲＷ等の光ディスク、ＤＶＤ等の相変化型光ディスク等で構成されてもよい。また、プログラムコードを、ネットワークを介してダウンロードしてもよい。
【００４１】
また、コンピュータが読み出したプログラムコードを実行することにより、上記実施の形態の機能が実現されるだけでなく、そのプログラムコードの指示に基づき、コンピュータ上で稼動しているＯＳ（オペレーティングシステム）等が実際の処理の一部または全部を行い、その処理によって前述した実施形態の機能が実現される場合も含まれることはいうまでもない。
【００４２】
さらに、記憶媒体から読み出されたプログラムコードが、コンピュータに挿入された機能拡張ボードやコンピュータに接続された機能拡張ユニットに備わるメモリに書き込まれた後、そのプログラムコードの指示に基づき、その機能拡張ボードや機能拡張ユニットに備わるＣＰＵ等が実際の処理の一部または全部を行い、その処理によって前述した実施形態の機能が実現される場合も含まれる。
【図面の簡単な説明】
【００４３】
【図１】本実施の形態に係る挙動シミュレーション方法に用いられる数値計算対象領域を概略的に示す図である。
【図２】図１の数値計算対象領域中の粒子が移動する状態を概略的に示す図である。
【図３】図１における粒子の運動状態計算処理の手順を示すフローチャートである。
【図４】本発明の実施の形態に係る計算装置としてのコンピュータの構成を概略的に示すブロック図である。
【符号の説明】
【００４４】
１ 数値計算対象領域
２ 空間要素
４ 流体
５ 粒子
【特許請求の範囲】
【請求項１】
容器内の流体中に分散する粒子の挙動シミュレーションを行う計算装置であって、
前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割手段と、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出手段と、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定手段と、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算手段と、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出手段と、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出手段と、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出手段とを備えることを特徴とする計算装置。
【請求項２】
前記設定手段は、前記粒子に対応する物理量を前記粒子に対応する前記空間要素に設定し、前記流体に対応する物理量を前記流体に対応する空間要素に設定することを特徴とする請求項１記載の計算装置。
【請求項３】
前記算出された前記粒子が受ける力は、前記数値計算対象領域の加速度と前記粒子の質量との積に依存した力を備えることを特徴とする請求項１又は２記載の計算装置。
【請求項４】
前記算出された前記粒子が受ける力は、前記粒子の電荷及び前記数値計算対象領域の電場の積に依存した力を備えることを特徴とする請求項１乃至３のいずれか１項に記載の計算装置。
【請求項５】
前記流体は、その内部に他の粒子を分散し、
前記算出された前記粒子が受ける力は、前記粒子と前記他の粒子との距離と相対速度に依存した力を備えることを特徴とする請求項１乃至４のいずれか１項に記載の計算装置。
【請求項６】
前記数値計算対象領域は、その中に固体領域が設定し、
前記算出された前記粒子が受ける力は、前記粒子と前記設定された固体領域との距離と前記粒子の速度に依存した力を備えることを特徴とする請求項１乃至５のいずれか１項に記載の計算装置。
【請求項７】
容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法であって、
前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割ステップと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出ステップと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定ステップと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算ステップと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出ステップと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出ステップと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出ステップとを備えることを特徴とする挙動シミュレーション方法。
【請求項８】
容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法をコンピュータにより実行させるプログラムであって
前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割モジュールと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出モジュールと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定モジュールと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算モジュールと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出モジュールと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出モジュールと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出モジュールとを備えることを特徴とするプログラム。
【請求項１】
容器内の流体中に分散する粒子の挙動シミュレーションを行う計算装置であって、
前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割手段と、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出手段と、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定手段と、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算手段と、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出手段と、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出手段と、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出手段とを備えることを特徴とする計算装置。
【請求項２】
前記設定手段は、前記粒子に対応する物理量を前記粒子に対応する前記空間要素に設定し、前記流体に対応する物理量を前記流体に対応する空間要素に設定することを特徴とする請求項１記載の計算装置。
【請求項３】
前記算出された前記粒子が受ける力は、前記数値計算対象領域の加速度と前記粒子の質量との積に依存した力を備えることを特徴とする請求項１又は２記載の計算装置。
【請求項４】
前記算出された前記粒子が受ける力は、前記粒子の電荷及び前記数値計算対象領域の電場の積に依存した力を備えることを特徴とする請求項１乃至３のいずれか１項に記載の計算装置。
【請求項５】
前記流体は、その内部に他の粒子を分散し、
前記算出された前記粒子が受ける力は、前記粒子と前記他の粒子との距離と相対速度に依存した力を備えることを特徴とする請求項１乃至４のいずれか１項に記載の計算装置。
【請求項６】
前記数値計算対象領域は、その中に固体領域が設定し、
前記算出された前記粒子が受ける力は、前記粒子と前記設定された固体領域との距離と前記粒子の速度に依存した力を備えることを特徴とする請求項１乃至５のいずれか１項に記載の計算装置。
【請求項７】
容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法であって、
前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割ステップと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出ステップと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定ステップと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算ステップと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出ステップと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出ステップと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出ステップとを備えることを特徴とする挙動シミュレーション方法。
【請求項８】
容器内の流体中に分散する粒子の挙動シミュレーションを行う挙動シミュレーション方法をコンピュータにより実行させるプログラムであって
前記容器内の領域を数値計算対象領域とし、当該数値計算対象領域を前記粒子よりも小さい空間要素に分割する領域分割モジュールと、現在時刻における前記粒子の位置・速度に基づいて、前記粒子が受ける力を算出する力算出モジュールと、前記算出された前記粒子が受ける力を前記粒子に対応する前記空間要素に設定する設定モジュールと、前記粒子及び前記流体を一体の流体として前記現在時刻から微小時間の間の流体計算を行う流体計算モジュールと、前記現在時刻から前記微小時間経過後の前記一体の流体の速度分布を算出する速度分布算出モジュールと、前記算出された速度分布のうち、前記粒子に対応する空間要素の速度を加重平均して前記現在時刻から前記微小時間経過後の前記粒子の速度を算出する粒子速度算出モジュールと、前記算出された粒子速度から前記微小時間経過後の前記粒子の位置を算出する粒子位置算出モジュールとを備えることを特徴とするプログラム。
【図１】
【図２】
【図３】
【図４】
【図２】
【図３】
【図４】
【公開番号】特開２００７−２８６９５５（Ｐ２００７−２８６９５５Ａ）
【公開日】平成１９年１１月１日（２００７．１１．１）
【国際特許分類】
【出願番号】特願２００６−１１４４５２（Ｐ２００６−１１４４５２）
【出願日】平成１８年４月１８日（２００６．４．１８）
【出願人】（０００００１００７）キヤノン株式会社 (59,756)
【Ｆターム（参考）】
【公開日】平成１９年１１月１日（２００７．１１．１）
【国際特許分類】
【出願日】平成１８年４月１８日（２００６．４．１８）
【出願人】（０００００１００７）キヤノン株式会社 (59,756)
【Ｆターム（参考）】
[ Back to top ]