WO2014045594A1 - 解析装置 - Google Patents
解析装置 Download PDFInfo
- Publication number
- WO2014045594A1 WO2014045594A1 PCT/JP2013/005581 JP2013005581W WO2014045594A1 WO 2014045594 A1 WO2014045594 A1 WO 2014045594A1 JP 2013005581 W JP2013005581 W JP 2013005581W WO 2014045594 A1 WO2014045594 A1 WO 2014045594A1
- Authority
- WO
- WIPO (PCT)
- Prior art keywords
- particle
- temperature
- particles
- unit
- force
- Prior art date
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume, or surface-area of porous materials
- G01N15/10—Investigating individual particles
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16C—COMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
- G16C10/00—Computational theoretical chemistry, i.e. ICT specially adapted for theoretical aspects of quantum chemistry, molecular mechanics, molecular dynamics or the like
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16C—COMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
- G16C99/00—Subject matter not provided for in other groups of this subclass
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume, or surface-area of porous materials
- G01N2015/0003—Determining electric mobility, velocity profile, average speed or velocity of a plurality of particles
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/22—Moulding
Landscapes
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Physics & Mathematics (AREA)
- Computing Systems (AREA)
- General Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Health & Medical Sciences (AREA)
- General Engineering & Computer Science (AREA)
- General Health & Medical Sciences (AREA)
- Mathematical Physics (AREA)
- Chemical & Material Sciences (AREA)
- Data Mining & Analysis (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Computational Mathematics (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- Dispersion Chemistry (AREA)
- Databases & Information Systems (AREA)
- Analytical Chemistry (AREA)
- Biochemistry (AREA)
- Algebra (AREA)
- Immunology (AREA)
- Pathology (AREA)
- Software Systems (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
解析装置100は、仮想空間に定義される粒子系の粒子の運動を支配する支配方程式を数値的に演算することによりその粒子系を解析する解析装置であって、粒子系の粒子のパラメータのひとつである温度を演算する温度演算部110と、粒子をその粒子について温度演算部110によって演算された温度の熱浴に浸すと仮定した場合の、その粒子に働く力を演算する力演算部122と、力演算部122によって演算された力に基づき粒子の状態を更新する状態更新部126と、を備える。
Description
本発明は、粒子系を解析する解析装置に関する。
古典力学や量子力学等を基に計算機を用いて物質科学全般の現象を探るための方法として、分子動力学法(Molecular Dynamics Method、以下MD法と称す)や、量子分子動力学法(Quantum Molecular Dynamics Method)や、MD法をマクロスケールの系を扱えるように発展させた繰り込み群分子動力学法(Renormalized Molecular Dynamics、以下RMD法と称す)に基づくシミュレーションが知られている。
解析対象の熱伝導現象について、MD法やRMD法では通常、格子振動(フォノン)による熱伝導しか扱えない。したがって自由電子の寄与が大きい金属などの熱伝導については、MD法やRMD法では実際とは乖離した解析結果が出ることが多い。
そこで本出願人は特許文献1および特許文献2において、いくつかの解決法を提案している。
John C. Tully、「Dynamics of gas-surface interaction: 3D generalized Langevin model applied to fcc and bcc surface」、J. Chem. Phys.、1975、73
特許文献1に記載の手法では、粒子に温度のパラメータを持たせ、有限体積法(Finite Volume Method、FVM)により熱伝導方程式を解くことで温度場を求めている。しかしながら、その温度場は粒子の速度の分散を反映していない場合が多く、適用範囲は比較的限られている。
特許文献2には、粒子の運動エネルギを確率変数に換算し、その勾配から、フーリエ則に従ってエネルギを受け渡す技術が開示されている。しかしながら、温度の勾配によってはエネルギが大きくなり、粒子の持つエネルギがマイナスになる場合がある。それを防ぐには時間刻みを小さく取ればよいのだが、そうすると計算負荷が増大しうる。
本発明はこうした課題に鑑みてなされたものであり、その目的は、複数の粒子を含む粒子系を使用したシミュレーションにおいて、熱伝導現象をより正確に扱えるようにする解析技術の提供にある。
本発明のある態様は、解析装置に関する。この解析装置は、仮想空間に定義される粒子系の粒子の運動を支配する支配方程式を数値的に演算することによりその粒子系を解析する解析装置であって、粒子系の粒子のパラメータのひとつである温度を演算する温度演算部と、粒子をその粒子について温度演算部によって演算された温度の熱浴に浸すと仮定した場合の、その粒子に働く力を演算する力演算部と、力演算部によって演算された力に基づき粒子の状態を更新する状態更新部と、を備える。
この態様によると、粒子系の熱伝導をより好適に扱うことができる。
なお、以上の構成要素の任意の組み合わせや、本発明の構成要素や表現を装置、方法、システム、コンピュータプログラム、コンピュータプログラムを格納した記録媒体などの間で相互に置換したものもまた、本発明の態様として有効である。
本発明によれば、複数の粒子を含む粒子系を使用したシミュレーションにおいて、熱伝導現象をより正確に扱うことができる。
以下、各図面に示される同一または同等の構成要素、部材、処理には、同一の符号を付するものとし、適宜重複した説明は省略する。
実施の形態に係る解析装置は、解析対象を複数の粒子を含む粒子系で記述し、粒子の運動方程式を数値的に演算することにより粒子系を解析する。その際、解析装置は熱伝導方程式を解いて温度場を求める。解析装置は、粒子の運動方程式に、求められた温度と実質的に等しい温度の熱浴に粒子が浸されると仮定したときに得られる粒子に働く力を導入する。具体的には、熱浴の粒子からのランダム力と、粘性力によるダンパーと、を与えるLangevin法(非特許文献1参照)が適用される。これにより、自由電子からの熱伝導への寄与が無視できない金属などの解析対象の熱伝導を再現する際の物理的な整合性を高めることができる。
(第1の実施の形態)
図1は、第1の実施の形態に係る解析装置100の機能および構成を示すブロック図である。ここに示す各ブロックは、ハードウエア的には、コンピュータのCPU(central processing unit)をはじめとする素子や機械装置で実現でき、ソフトウエア的にはコンピュータプログラム等によって実現されるが、ここでは、それらの連携によって実現される機能ブロックを描いている。したがって、これらの機能ブロックはハードウエア、ソフトウエアの組合せによっていろいろなかたちで実現できることは、本明細書に触れた当業者には理解されるところである。
図1は、第1の実施の形態に係る解析装置100の機能および構成を示すブロック図である。ここに示す各ブロックは、ハードウエア的には、コンピュータのCPU(central processing unit)をはじめとする素子や機械装置で実現でき、ソフトウエア的にはコンピュータプログラム等によって実現されるが、ここでは、それらの連携によって実現される機能ブロックを描いている。したがって、これらの機能ブロックはハードウエア、ソフトウエアの組合せによっていろいろなかたちで実現できることは、本明細書に触れた当業者には理解されるところである。
本実施の形態ではMD法に倣って粒子系を解析する場合について説明するが、RMD法に倣って粒子系を解析する場合や、DEM(Distinct Element Method)やSPH(Smoothed Particle Hydrodynamics)やMPS(Moving Particle Semi-implicit)などの粒子法に倣って粒子系を解析する場合にも、本実施の形態に係る技術的思想を適用できることは本明細書に触れた当業者には明らかである。
解析装置100は入力装置102およびディスプレイ104と接続される。入力装置102は、解析装置100で実行される処理に関係するユーザの入力を受けるためのキーボード、マウスなどであってもよい。入力装置102は、インターネットなどのネットワークやCD、DVDなどの記録媒体から入力を受けるよう構成されていてもよい。
解析装置100は、粒子系取得部108と、温度付与部134と、繰り返し演算部120と、表示制御部118と、粒子データ保持部114と、を備える。
粒子系取得部108は、入力装置102を介してユーザから取得する入力情報に基づき、1、2または3次元の仮想空間内に定義されるN(Nは自然数)個の粒子からなる粒子系のデータを取得する。粒子系の粒子は分子または原子に対応してもよい。
粒子系取得部108は、入力情報に基づき仮想空間内にN個の粒子を配置し、配置されたそれぞれの粒子に速度を付与する。すなわち、粒子系取得部108は、粒子系に初期位置および初期速度を付与する。粒子系取得部108は、配置された粒子を特定する粒子IDと、その粒子の位置と、その粒子の速度と、を対応付けて粒子データ保持部114に登録する。
温度付与部134は、入力装置102を介してユーザから取得する入力情報に基づき、粒子系取得部108によって取得された粒子系の粒子に温度を付与する。このように付与される温度は粒子のパラメータのひとつである。例えば、温度付与部134は、ディスプレイ104を介してユーザに、粒子系の粒子の温度の初期値の入力を求める。温度付与部134は、入力された温度の初期値を粒子IDと対応付けて粒子データ保持部114に登録する。
以下では粒子系の粒子は全て同質または同等なものとして設定され、かつ、ポテンシャルエネルギ関数は2体のポテンシャルであって粒子によらずに同じ形を有するものとして設定される場合について説明する。しかしながら、他の場合にも本実施の形態に係る技術的思想を適用できることは、本明細書に触れた当業者には明らかである。
繰り返し演算部120は、粒子データ保持部114によって保持されるデータが表す粒子系の各粒子の運動を支配する支配方程式を数値的に演算する。特に繰り返し演算部120は、離散化された粒子の運動方程式にしたがった繰り返し演算を行う。
繰り返し演算部120は、温度演算部110と、力演算部122と、粒子状態演算部124と、状態更新部126と、終了条件判定部128と、を含む。
繰り返し演算部120は、温度演算部110と、力演算部122と、粒子状態演算部124と、状態更新部126と、終了条件判定部128と、を含む。
温度演算部110は、連続体近似を用いて粒子系の各粒子の温度を演算する。特に温度演算部110は、離散化された熱伝導方程式に基づき粒子の温度を演算する。温度演算部110は、ボロノイ分割部112と、熱伝導演算部116と、を有する。
ボロノイ分割部112は、粒子系が定義されている仮想空間の部分をボロノイ分割する。すなわち、ボロノイ分割部112は、各粒子の位置に基づいて仮想空間内にボロノイ多面体を生成する。ボロノイ分割部112はまず、各粒子を頂点として仮想空間の部分を3次元デローニ(ドロネー)分割する。次にボロノイ分割部112は、デローニ分割の結果得られる四面体要素からボロノイ多面体要素を生成する。これにより、仮想空間の部分は粒子の位置を母点としてボロノイ分割される。
熱伝導演算部116は、ボロノイ分割部112によって生成されたボロノイ多面体を単位として離散化された熱伝導方程式に基づいて各粒子の温度を演算する。熱伝導演算部116では、有限体積法を用いた温度場の解析が実行されてもよい。熱伝導演算部116は、温度場の解析において、ボロノイ分割部112によって生成されたボロノイ多面体の面積や体積、および所定の微小な時間刻みΔtだけ前の時刻の温度場の情報(すなわち、繰り返し演算における1サイクル前の温度場の情報)を使用する。特に熱伝導演算部116は、ボロノイ多面体のそれぞれをFVMのひとつのコントロールボリュームとして解析する。
以下、熱伝導演算部116で使用される離散化された熱伝導方程式の導出過程を示す。
熱伝導方程式は以下の微分方程式(式1)で与えられる。
ここで、ρは密度、Cvは比熱、Tは温度、tは時間、λは熱伝導率、Qは単位体積当たりの発熱量である。
熱伝導方程式は以下の微分方程式(式1)で与えられる。
式3をボロノイ多面体を単位として離散化すると、
となる。ここで、ΔViはi番目の粒子の位置を母点とするボロノイ多面体の体積、ΔSijはi番目の粒子とj番目の粒子との間にあるボロノイ面の面積、rijはi番目の粒子とj番目の粒子との距離である。また、Ti
nは時間に関する繰り返し演算のn回目におけるi番目の粒子の温度(i番目の粒子が属するボロノイ多面体の温度とも言える)、λijはi番目の粒子とj番目の粒子との間の熱伝導率である。式4は、熱伝導演算部116で使用される離散化された熱伝導方程式である。式4により、ボロノイ多面体の体積および面積と、n回目の演算により決定された各粒子の温度と、粒子間距離と、から、n+1回目の演算におけるi番目の粒子の温度を決定することができる。なお、Ti
0には、温度付与部134によって付与された初期温度を使用する。
力演算部122は、粒子をその粒子について温度演算部110によって演算された温度の熱浴に浸すと仮定した場合の、その粒子に働く力を演算する。力演算部122は、粒子間作用演算部130と、熱浴作用演算部132と、を有する。
粒子間作用演算部130は、粒子データ保持部114によって保持される粒子系のデータを参照し、粒子系の各粒子について、粒子間の距離に基づきその粒子に働く力を演算する。粒子間作用演算部130は、粒子系のi番目(1≦i≦N)の粒子について、そのi番目の粒子との距離が所定のカットオフ距離よりも小さな粒子(以下、近接粒子と称す)を決定する。
粒子間作用演算部130は、各近接粒子について、その近接粒子とi番目の粒子との間のポテンシャルエネルギ関数およびその近接粒子とi番目の粒子との距離に基づいて、その近接粒子がi番目の粒子に及ぼす力を演算する。特に粒子間作用演算部130は、その近接粒子とi番目の粒子との距離の値におけるポテンシャルエネルギ関数のグラジエント(Gradient)の値から力を算出する。粒子間作用演算部130は、近接粒子がi番目の粒子に及ぼす力を全ての近接粒子について足し合わせることによって、i番目の粒子に働く力を算出する。粒子間作用演算部130によって算出される力は、粒子間の相互作用に基づく力である。
i番目の粒子について温度演算部110によって演算された温度Tiで実質的に定温の熱浴にi番目の粒子を浸すと仮定した場合、Langevin法によるとi番目の粒子に以下の2つの力が働く。
(1)ダンパーによる力(粘性力)
熱浴の粘性が粒子に作用する力である。粘性力の減衰定数αは、デバイ周波数ωDと粒子の質量mとを用いて、
と表される。以下の表は、代表的な金属物質の減衰定数αの値を示す。この表では粒子を原子に対応付けている。
この表に示される通り、粒子が金属原子に対応する場合、減衰定数αは1.0x10-12(kg/s)のオーダーとなる。なお、デバイ周波数ωDは粒子の質量に依存するので、粒子の質量が原子の質量のβ倍になると、減衰定数αはβ0.5倍となる。例えば、鉄の物性を持ちながら質量が鉄原子の100倍大きい粒子を用いる場合、表1に従うと減衰定数αは2.99x10-11(kg/s)となる。
熱浴の粘性が粒子に作用する力である。粘性力の減衰定数αは、デバイ周波数ωDと粒子の質量mとを用いて、
熱浴作用演算部132は、粒子系のi番目(1≦i≦N)の粒子について、式5および式6に基づき粘性力とランダム力とを演算する。熱浴作用演算部132は、粒子間の相互作用に基づくi番目の粒子に働く力に、i番目の粒子について演算された粘性力およびランダム力を加えることによって、i番目の粒子に働くトータルの力を演算する。i番目の粒子に働くトータルの力Fiは以下の式7で表される。
ここで、φijはi番目の粒子とj番目の粒子との間のポテンシャルエネルギ関数、viはi番目の粒子の速度、Frandomは標準偏差σを有するランダム力である。記号の上の矢印はベクトル量であることを示す。
粒子状態演算部124は粒子データ保持部114に保持される粒子系のデータを参照し、粒子系の各粒子について、離散化された粒子の運動方程式に熱浴作用演算部132によって演算されたトータルの力を適用することによって粒子の位置および速度のうちの少なくともひとつを演算する。本実施の形態では、粒子状態演算部124は粒子の位置および速度の両方を演算する。
粒子状態演算部124は、熱浴作用演算部132によって演算されたトータルの力を含む離散化された粒子の運動方程式から粒子の速度を演算する。粒子状態演算部124は、粒子系のi番目の粒子について、蛙跳び法やオイラー法などの所定の数値解析の手法に基づき時間刻みΔtを使用して離散化された粒子の運動方程式に、i番目の粒子について熱浴作用演算部132によって演算されたトータルの力を代入することによって、粒子の速度を演算する。この演算には以前の繰り返し演算のサイクルで演算された粒子の速度が使用される。
粒子状態演算部124は、演算された粒子の速度に基づいて粒子の位置を算出する。粒子状態演算部124は、粒子系のi番目の粒子について、所定の数値解析の手法に基づき時間刻みΔtを使用して離散化された粒子の位置と速度の関係式に、演算された粒子の速度を適用することによって、粒子の位置を演算する。この演算には以前の繰り返し演算のサイクルで演算された粒子の位置が使用される。
状態更新部126は、粒子状態演算部124における演算結果に基づき粒子系の各粒子の状態を更新する。状態更新部126は、粒子データ保持部114に保持される粒子系の各粒子の位置および速度のそれぞれを、粒子状態演算部124によって演算された位置および速度で更新する。
終了条件判定部128は、繰り返し演算部120における繰り返し演算を終了すべきか否かを判定する。繰り返し演算を終了すべき終了条件は、例えば繰り返し演算が所定の回数行われたことや、外部から終了の指示を受け付けたことや、粒子系が定常状態に達したことである。終了条件判定部128は、終了条件が満たされる場合、繰り返し演算部120における繰り返し演算を終了させる。終了条件判定部128は、終了条件が満たされない場合、処理を温度演算部110に戻す。すると温度演算部110は、状態更新部126によって更新された粒子の位置で再び温度を演算する。
表示制御部118は、粒子データ保持部114に保持されるデータが表す粒子系の各粒子の位置、速度、温度に基づき、ディスプレイ104に粒子系の時間発展の様子やある時刻における状態を表示させる。この表示は、静止画または動画の形式で行われてもよい。
図2は、粒子データ保持部114の一例を示すデータ構造図である。粒子データ保持部114は、粒子IDと、粒子の位置と、粒子の速度と、粒子の温度と、を対応付けて保持する。
上述の実施の形態において、保持部の例は、ハードディスクやメモリである。また、本明細書の記載に基づき、各部を、図示しないCPUや、インストールされたアプリケーションプログラムのモジュールや、システムプログラムのモジュールや、ハードディスクから読み出したデータの内容を一時的に記憶するメモリなどにより実現できることは本明細書に触れた当業者には理解されるところである。
以上の構成による解析装置100の動作を説明する。
図3は、解析装置100における一連の処理の一例を示すフローチャートである。解析装置100は、粒子系の初期状態すなわち各粒子の初期位置、初期速度および初期温度を決定する(S202)。解析装置100は、各粒子の位置に基づきボロノイ解析を実行し、ボロノイ多面体を生成する(S204)。解析装置100は、FVMを使用して温度場を解析し(S206)、各粒子の温度を更新する。解析装置100は、粒子間のポテンシャルエネルギ関数に基づく各粒子に働く力を演算する(S208)。解析装置100は、ステップS208で演算された力に、粘性力およびランダム力を付加する(S210)。解析装置100は、ステップS210で演算された力を含む粒子の運動方程式から粒子の速度と位置を演算する(S212)。解析装置100は、粒子データ保持部114に保持される粒子の位置、速度を演算された位置、速度で更新する(S214)。解析装置100は、終了条件が満たされるか否かを判定する(S216)。終了条件が満たされない場合(S216のN)、処理はステップS204に戻される。終了条件が満たされる場合(S216のY)、処理は終了する。
図3は、解析装置100における一連の処理の一例を示すフローチャートである。解析装置100は、粒子系の初期状態すなわち各粒子の初期位置、初期速度および初期温度を決定する(S202)。解析装置100は、各粒子の位置に基づきボロノイ解析を実行し、ボロノイ多面体を生成する(S204)。解析装置100は、FVMを使用して温度場を解析し(S206)、各粒子の温度を更新する。解析装置100は、粒子間のポテンシャルエネルギ関数に基づく各粒子に働く力を演算する(S208)。解析装置100は、ステップS208で演算された力に、粘性力およびランダム力を付加する(S210)。解析装置100は、ステップS210で演算された力を含む粒子の運動方程式から粒子の速度と位置を演算する(S212)。解析装置100は、粒子データ保持部114に保持される粒子の位置、速度を演算された位置、速度で更新する(S214)。解析装置100は、終了条件が満たされるか否かを判定する(S216)。終了条件が満たされない場合(S216のN)、処理はステップS204に戻される。終了条件が満たされる場合(S216のY)、処理は終了する。
温度演算部110では連続体近似により各粒子の温度が導出されるので、温度演算部110によって演算される粒子の温度は、本来の温度の定義である粒子の速度の分散と大きく異なることが多い。本発明者はそのような非整合性を軽減または除去するために、温度演算部110によって温度を求めた後、その温度に由来する運動エネルギを粒子の運動に反映させることに想到した。ここで、例えば温度スケーリングなどで粒子の速度を強制的に温度に対応する速度に変更することが考えられるが、これは運動を強制的に拘束するため、非物理的である。
そこで、本実施の形態に係る解析装置100では、粒子を温度演算部110によって演算された温度の熱浴に浸すと仮定することで、運動方程式の力の項を温度に基づき修正している。これにより、温度演算部110によって演算された温度を粒子の速度場に反映させることができ、温度演算部110によって演算された温度場により自然な形で導くことができる。その結果、物理的な矛盾がより少ないモデルを提供できる。
本発明者は、本実施の形態に係る手法の確認計算を行った。本実施の形態が使用するMD法では、基本的に熱伝導として粒子の格子振動のみしか扱えないため、自由電子の寄与が反映されない。したがって、MD法で金属の解析対象を扱う場合、すなわち粒子系の粒子が金属の粒子を模するように粒子についての物質定数(例えば、デバイ温度やデバイ周波数や原子量、および、熱伝導方程式における密度や比熱や熱伝導率)が設定される場合、本実施の形態に係る手法が非常に有効である。
図4は、確認計算で使用された粒子系300を示す模式図である。粒子系300は金属の棒を模している。この棒の両端の温度は0(K)で固定され、初期の温度分布は以下の式8により与えられる。
ここで、Lは棒の長さ、rは端からの距離、T(r)は距離rにおける温度である。
図5は、本実施の形態に係る手法を使用した計算結果を示すグラフである。図6は、本実施の形態に係る手法を使用しない場合の計算結果を示すグラフである。本実施の形態に係る手法によると、粒子系300が時間発展を始めてから0.3(ns)経過後、0.6(ns)経過後、および0.9(ns)経過後のいずれにおいても、温度分布の計算値(黒点で表示)と理論値(実線で表示)とが良く一致していることが分かる。これに対して本実施の形態に係る手法を使用しない通常のMD法の場合、理論値と比較して熱の拡散が非常に遅くなっていることが分かる。
(第2の実施の形態)
第2の実施の形態では、粒子系の構造の変化、すなわち粒子配置の変化に伴う発熱を考慮した場合について説明する。第2の実施の形態に係る解析装置100は、図1と同様の構成を有する。以下、第1の実施の形態との相違点を中心に説明する。
第2の実施の形態では、粒子系の構造の変化、すなわち粒子配置の変化に伴う発熱を考慮した場合について説明する。第2の実施の形態に係る解析装置100は、図1と同様の構成を有する。以下、第1の実施の形態との相違点を中心に説明する。
熱伝導演算部116は、第1の実施の形態と同様に、ボロノイ分割部112によって生成されたボロノイ多面体を単位として離散化された熱伝導方程式に基づいて各粒子の温度を演算する。熱伝導演算部116で使用される離散化された熱伝導方程式は、第1の実施の形態と同様に、式4で表される。本実施の形態では、粒子系の構造の変化に伴う発熱を考慮するため、式4の発熱量Qは以下の式11で表される。
ここで、φは粒子間のポテンシャルエネルギ関数、Kiはi番目の粒子の運動エネルギ、nは時間に関する繰り返し演算の回数、Δtは時間刻み、Frijはi番目の粒子とj番目の粒子との間で生じる摩擦力、vrijはi番目の粒子とj番目の粒子との相対速度である。式11は、粒子系の変形により粒子のもつ全エネルギが増加すれば発熱することを示している。式4および式11により、n+1回目の演算におけるi番目の粒子の温度を決定することができる。
本実施の形態に係る解析装置100によれば、第1の実施の形態と同様に、温度演算部110によって演算された温度を粒子の速度場に反映させることができ、温度演算部110によって演算された温度場により自然な形で導くことができる。その結果、物理的な矛盾がより少ないモデルを提供できる。加えて、本実施の形態に係る解析装置100によれば、温度演算部110(熱伝導演算部116)により演算される温度に、粒子系の構造の変化による発熱を反映させることができる。そのため、粒子系の構造に変形が生じる場合についても、物理的な矛盾がより少ないモデルを提供できる。これにより、例えば塑性加工など金属変形に伴う発熱が関与する現象をより正確にシミュレートすることができ、加工時の温度上昇を予測することができる。
本発明者は、本実施の形態に係る手法の確認計算を行った。図7(a)~(f)は、本実施の形態に係る手法を使用した計算結果を示す模式図である。粒子系400は、0.6mm(X方向)×0.6mm(Y方向)×0.95mm(Z方向)の金属の角材を模している。図7(a)~(f)は、この粒子系400の上下1層にそれぞれ応力50GPa相当の引っ張り力を加えた場合の時間変化を示している。図7(a)は、粒子系400に引っ張り力を加えた瞬間を示している。図7(b)、(c)、(d)、(e)、(f)はそれぞれ、粒子系400が時間発展を始めてから25μs経過後、50μs経過後、75μs経過後、100μs経過後、125μs経過後の状態を示す。この結果から、引っ張り力が加えられ、粒子系400の構造が変化したことにより、粒子系400の温度が上昇していることが分かる。これは、粒子系400の構造が変化すると熱が発生するという知見と合致する。
以上、実施の形態に係る解析装置100の構成と動作について説明した。この実施の形態は例示であり、その各構成要素や各処理の組み合わせにいろいろな変形例が可能なこと、またそうした変形例も本発明の範囲にあることは当業者に理解されるところである。
実施の形態では、繰り返し演算部120において粒子の位置と速度の両方を演算する場合について説明したが、これに限られない。例えば、数値解析の手法にはVerlet法のように、粒子の位置を演算する際に粒子に働く力から粒子の位置を直接演算し、粒子の速度は陽に計算しなくてもよい手法もあり、本実施の形態に係る技術的思想をそのような手法に適用してもよい。
100 解析装置、 102 入力装置、 104 ディスプレイ、 108 粒子系取得部、 114 粒子データ保持部、 118 表示制御部、 120 繰り返し演算部。
本発明によれば、複数の粒子を含む粒子系を使用したシミュレーションにおいて、熱伝導現象をより正確に扱うことができる。
Claims (5)
- 仮想空間に定義される粒子系の粒子の運動を支配する支配方程式を数値的に演算することによりその粒子系を解析する解析装置であって、
粒子系の粒子のパラメータのひとつである温度を演算する温度演算部と、
粒子をその粒子について前記温度演算部によって演算された温度の熱浴に浸すと仮定した場合の、その粒子に働く力を演算する力演算部と、
前記力演算部によって演算された力に基づき粒子の状態を更新する状態更新部と、を備えることを特徴とする解析装置。 - 前記温度演算部は、
粒子の位置に基づいて仮想空間にボロノイ多面体を生成するボロノイ処理部と、
前記ボロノイ処理部によって生成されたボロノイ多面体を単位として離散化された熱伝導方程式に基づいて粒子の温度を演算する熱伝導演算部と、を含むことを特徴とする請求項1に記載の解析装置。 - 前記熱伝導方程式は、粒子系の変形による発熱を表す項を含むことを特徴とする請求項2に記載の解析装置。
- 粒子系の粒子が金属の粒子を模するように粒子についての物質定数が設定されることを特徴とする請求項1から3のいずれかに記載の解析装置。
- 仮想空間に定義される粒子系の粒子の運動を支配する支配方程式を数値的に演算することによりその粒子系を解析する機能をコンピュータに実現させるコンピュータプログラムであって、
粒子系の粒子のパラメータのひとつである温度を演算する機能と、
粒子をその粒子について演算された温度の熱浴に浸すと仮定した場合の、その粒子に働く力を演算する機能と、
演算された力に基づき粒子の状態を更新する機能と、を前記コンピュータに実現させることを特徴とするコンピュータプログラム。
Priority Applications (3)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2014536606A JP6129193B2 (ja) | 2012-09-21 | 2013-09-20 | 解析装置 |
EP13838490.4A EP2899655A4 (en) | 2012-09-21 | 2013-09-20 | ANALYSIS DEVICE |
US14/658,408 US20150186573A1 (en) | 2012-09-21 | 2015-03-16 | Analyzing device |
Applications Claiming Priority (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2012208720 | 2012-09-21 | ||
JP2012-208720 | 2012-09-21 |
Related Child Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
US14/658,408 Continuation US20150186573A1 (en) | 2012-09-21 | 2015-03-16 | Analyzing device |
Publications (1)
Publication Number | Publication Date |
---|---|
WO2014045594A1 true WO2014045594A1 (ja) | 2014-03-27 |
Family
ID=50340835
Family Applications (2)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
PCT/JP2013/003501 WO2014045492A1 (ja) | 2012-09-21 | 2013-06-04 | 解析装置 |
PCT/JP2013/005581 WO2014045594A1 (ja) | 2012-09-21 | 2013-09-20 | 解析装置 |
Family Applications Before (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
PCT/JP2013/003501 WO2014045492A1 (ja) | 2012-09-21 | 2013-06-04 | 解析装置 |
Country Status (4)
Country | Link |
---|---|
US (1) | US20150186573A1 (ja) |
EP (1) | EP2899655A4 (ja) |
JP (1) | JP6129193B2 (ja) |
WO (2) | WO2014045492A1 (ja) |
Families Citing this family (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US10311176B2 (en) | 2015-09-03 | 2019-06-04 | Sumitomo Heavy Industries, Ltd. | Simulation method, simulation apparatus, and simulation program |
JP6679161B2 (ja) * | 2015-09-03 | 2020-04-15 | 住友重機械工業株式会社 | シミュレーション方法、シミュレーション装置、及びシミュレーションプログラム |
JP6910729B2 (ja) * | 2017-10-06 | 2021-07-28 | 住友重機械工業株式会社 | シミュレーション方法、シミュレーション装置、及びプログラム |
CN107967405A (zh) * | 2017-11-27 | 2018-04-27 | 中国计量大学 | 提高分子动力学计算效率的方法 |
CN112069579B (zh) * | 2020-09-04 | 2022-08-05 | 华能澜沧江水电股份有限公司 | 基于dem数字地形分析的土石坝变形震害定量评估方法 |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPS61132205A (ja) * | 1984-12-03 | 1986-06-19 | Kawasaki Steel Corp | 珪素鋼板の冷間圧延方法 |
JP2001334305A (ja) * | 2000-05-22 | 2001-12-04 | Nkk Corp | 鋼板の仕上げ温度制御方法 |
JP2010170309A (ja) | 2009-01-22 | 2010-08-05 | Sumitomo Heavy Ind Ltd | シミュレーション方法及びプログラム |
JP2010230331A (ja) * | 2009-03-25 | 2010-10-14 | Neturen Co Ltd | 高周波焼入れシミュレーション装置 |
JP2011175327A (ja) * | 2010-02-23 | 2011-09-08 | Sumitomo Heavy Ind Ltd | シミュレーション方法及びプログラム |
JP2012128490A (ja) * | 2010-12-13 | 2012-07-05 | Sumitomo Heavy Ind Ltd | 解析装置および解析方法 |
JP2012150673A (ja) * | 2011-01-19 | 2012-08-09 | Sumitomo Heavy Ind Ltd | 解析装置およびシミュレーション方法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP5052985B2 (ja) * | 2007-07-31 | 2012-10-17 | 住友重機械工業株式会社 | 分子シミュレーション方法、分子シミュレーション装置、分子シミュレーションプログラム、及び該プログラムを記録した記録媒体 |
-
2013
- 2013-06-04 WO PCT/JP2013/003501 patent/WO2014045492A1/ja active Application Filing
- 2013-09-20 JP JP2014536606A patent/JP6129193B2/ja active Active
- 2013-09-20 EP EP13838490.4A patent/EP2899655A4/en not_active Withdrawn
- 2013-09-20 WO PCT/JP2013/005581 patent/WO2014045594A1/ja active Application Filing
-
2015
- 2015-03-16 US US14/658,408 patent/US20150186573A1/en not_active Abandoned
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPS61132205A (ja) * | 1984-12-03 | 1986-06-19 | Kawasaki Steel Corp | 珪素鋼板の冷間圧延方法 |
JP2001334305A (ja) * | 2000-05-22 | 2001-12-04 | Nkk Corp | 鋼板の仕上げ温度制御方法 |
JP2010170309A (ja) | 2009-01-22 | 2010-08-05 | Sumitomo Heavy Ind Ltd | シミュレーション方法及びプログラム |
JP2010230331A (ja) * | 2009-03-25 | 2010-10-14 | Neturen Co Ltd | 高周波焼入れシミュレーション装置 |
JP2011175327A (ja) * | 2010-02-23 | 2011-09-08 | Sumitomo Heavy Ind Ltd | シミュレーション方法及びプログラム |
JP2012128490A (ja) * | 2010-12-13 | 2012-07-05 | Sumitomo Heavy Ind Ltd | 解析装置および解析方法 |
JP2012150673A (ja) * | 2011-01-19 | 2012-08-09 | Sumitomo Heavy Ind Ltd | 解析装置およびシミュレーション方法 |
Non-Patent Citations (2)
Title |
---|
JOHN C. TULLY: "Dynamics of gas-surface interaction: 3D generalized Langevin model applied to fcc and bcc surface", J. CHEM. PHYS., 1975, pages 73 |
See also references of EP2899655A4 |
Also Published As
Publication number | Publication date |
---|---|
EP2899655A4 (en) | 2016-04-27 |
JPWO2014045594A1 (ja) | 2016-08-18 |
EP2899655A1 (en) | 2015-07-29 |
JP6129193B2 (ja) | 2017-05-17 |
WO2014045492A1 (ja) | 2014-03-27 |
US20150186573A1 (en) | 2015-07-02 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Kang et al. | A direct-forcing immersed boundary method for the thermal lattice Boltzmann method | |
CN107016154A (zh) | 在物理坐标中有效地求解具有模态阻尼的情况下的结构动力学问题 | |
Servin et al. | Examining the smooth and nonsmooth discrete element approaches to granular matter | |
JP6129193B2 (ja) | 解析装置 | |
Peer et al. | Prescribed velocity gradients for highly viscous SPH fluids with vorticity diffusion | |
US20150356217A1 (en) | Lattice Boltzmann Collision Operators Enforcing Isotropy and Galilean Invariance | |
Venkatesan et al. | A three-field local projection stabilized formulation for computations of Oldroyd-B viscoelastic fluid flows | |
Jun et al. | Assessment of the cubic Fokker–Planck–DSMC hybrid method for hypersonic rarefied flows past a cylinder | |
Mollon | A multibody meshfree strategy for the simulation of highly deformable granular materials | |
JP6679161B2 (ja) | シミュレーション方法、シミュレーション装置、及びシミュレーションプログラム | |
Lorenz et al. | Combined lattice–Boltzmann and rigid-body method for simulations of shear-thickening dense suspensions of hard particles | |
EP2808812A1 (en) | Analysis device and simulation method | |
WO2014045416A1 (ja) | シミュレーションプログラム、シミュレーション方法及びシミュレーション装置 | |
Cinquegrana et al. | Validation of a new fluid—structure interaction framework for non-linear instabilities of 3D aerodynamic configurations | |
JP6053418B2 (ja) | 解析方法および解析装置 | |
Jarosch | Icetools: A full Stokes finite element model for glaciers | |
JP5930987B2 (ja) | 解析装置およびコンピュータプログラム | |
Kononenko et al. | Machine learning and finite element method for physical systems modeling | |
Kalarakis et al. | Lattice‐Boltzmann and meshless point collocation solvers for fluid flow and conjugate heat transfer | |
Hashemi et al. | Comparative study of momentum-exchange and smoothed profile methods in Lattice Boltzmann method | |
JP6091402B2 (ja) | 解析装置および解析方法 | |
Blom et al. | Partitioned Fluid–Structure–Acoustics Interaction on Distributed Data: Numerical Results and Visualization | |
Vasyliv et al. | Simulating incompressible flow on moving meshfree grids | |
Kang et al. | An accurate and efficient method for automatic deformation of unstructured polyhedral grids to simulate the flow induced by the motion of solid objects | |
JP6029457B2 (ja) | 解析装置および解析方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
121 | Ep: the epo has been informed by wipo that ep was designated in this application |
Ref document number: 13838490 Country of ref document: EP Kind code of ref document: A1 |
|
ENP | Entry into the national phase |
Ref document number: 2014536606 Country of ref document: JP Kind code of ref document: A |
|
WWE | Wipo information: entry into national phase |
Ref document number: 2013838490 Country of ref document: EP |
|
NENP | Non-entry into the national phase |
Ref country code: DE |