TECHNOLOGY EXPLAINER
分子動力学
― 原子をニュートンの運動方程式で動かすと、何が分かり、どこで限界が来るか
金属の中で原子はどう拡散するのか。高分子の鎖は温度で、どう動きやすくなるのか。固体電解質の中を、イオンはどう進むのか。分子動力学（MD）は、原子を小さな玉と見なし、各原子にかかる力からニュートンの運動方程式で「ほんの少しずつ」動かして、その動きを時間を追って観察する計算です。この記事では、図と小さな計算例を使って、MDの仕組み、力場（原子のあいだの力の決め方）の役割、時間と大きさの限界、何が分かるのかを、素材技術者の目線で整理します。
- 3行でわかる／読む前に／シリーズの中での位置
- まず身近な感覚：原子を動かす「スロー再生」
- 力はどこから来るか：エネルギーの坂
- 時間刻み：なぜ約1フェムト秒なのか
- 力場の役割と種類
- 温度・圧力を保つ方法と、小さな箱の工夫
- 時間と大きさの限界
- MDから何が分かるか
- 研究の進め方と、確かめるところ
- 素材技術者から見て
- 限界と、確認できなかったこと／まとめ
- 評価・標準の読み方
- 用語集
- 参考資料（一次情報）
- 主張と出典の対応表（監査）
事実＝公表資料に記載のある内容（出典リンクあり）
本稿の計算＝明示した前提をもとに本記事が算出した値
未確定／将来＝計画・目標・草案・改訂中・策定中・見通しなど、確定していない事項
このほか、素材技術者向けの読み取りは「解説」として区別しています。推測は太字で示します。
1. 3行でわかる／読む前に／シリーズの中での位置
- 分子動力学（MD）は、原子を玉と見なし、ニュートンの運動方程式（力＝質量×加速度）で、ごく短い時間ずつ動かしていく計算です。1ステップは、代表的なソフトウェア LAMMPS の既定値で、単位の選び方によって 1 fs（フェムト秒、10⁻¹⁵秒）または 0.001 ps（＝1 fs）です事実［参考資料9］
- 結果の良し悪しは、「原子のあいだの力をどう決めるか（力場）」でほぼ決まります。また、1 ms（ミリ秒）を動かすには、1 fs のステップを 10¹² 回くり返す計算が要ります本稿の計算。時間と大きさには、はっきりした限界があります
- MDで分かるのは、拡散のしやすさ、熱の伝わり方、変形の仕方、構造の変化など、原子の動きに由来する量です。ただし、実測と比べて初めて意味が出ます。AIが材料を発明するのではなく、MDは「原子の動きを計算で観察する道具」です（本記事による解説）
必要な前提は、「力が決まれば、物体の加速度が決まる」というニュートンの運動方程式だけです。それ以外の用語（ポテンシャル、力場、アンサンブルなど）は、初めて出るところで説明します。この記事は単独でも読めるよう書いています。
本シリーズの共通の考え方も確認しておきます。材料探索は、広大な地図の中から宝のありかを探す作業です。計算は「安く速く候補を出す」、実験は「本当にそうか確かめる」、AIは「探す場所を絞る」という分担です。MDは、この中の「計算」の一つで、原子の動きを時間を追って見られることが特徴です。
電子の振る舞いから性質を計算する方法は、本シリーズ「密度汎関数理論の基礎」編と「量子化学・材料シミュレーション」編で扱います。この記事（分子動力学）は、原子の動きを追う方法です。次に「機械学習原子間ポテンシャル」編で、力の計算を速くする方法を、「マルチスケールシミュレーションとデジタルツイン」編で、電子から部材までのつなぎ方を扱います。「不確かさと適用範囲」編と合わせて読むと、MDの結果をどこまで信じてよいかが分かります。
2. まず身近な感覚：原子を動かす「スロー再生」
ビリヤード台を、ものすごいスロー再生で撮った動画を想像してください。球が1個ずつ、どの向きに、どれだけ動くかが、コマ送りで分かります。MDは、原子でこれをやります。箱の中に原子を並べ、「いま、各原子にどんな力がかかっているか」を計算し、その力で原子をほんの少し動かし、また力を計算する。これをひたすらくり返すと、原子の動きの「動画」ができます。
この動画では、各原子の位置が時刻ごとに分かります。そこから、「どれだけ動き回ったか（拡散）」「原子の並びがどう変わったか」「箱を引き伸ばしたとき、どれだけの力が要ったか」といった、材料の性質に結びつく量を計算できます。
元のSVG・数値・出典を確認する
元SVGの説明：概念図（ベクター作図）。MDの1ステップ。①いまの位置から力を求め、②ニュートンの運動方程式（F＝m a）で加速度を出し、③ごく短い時間 Δt だけ速度と位置を更新し、④記録して①へ戻る。
ここで使う計算の手順を、積分法（運動方程式を少しずつ進める計算の決まり）といいます。もっとも標準的なのは、Verlet（1967）の論文で使われた方法に由来するベルレ法で、LAMMPS の既定の方法は、その速度版（速度ベルレ法）です事実［参考資料7・8］。MDの考え方そのものは古く、Alder と Wainwright（1957）が硬い球の系の固体と液体の転移を調べ、Rahman（1964）が液体アルゴンの原子の動きを計算し、McCammon、Gelin、Karplus（1977）がタンパク質（BPTI）の動きを計算しました事実［参考資料1・2・3］。
MDの原子は、電子を含めた全体を1つの玉に見立てたものです。電子の状態は、力の計算（力場）の中に「ひとまとめ」にして入っています。原子の動きは、量子力学ではなく、ニュートンの運動方程式（古典力学）で追います。この近似は、水素のように軽い原子や、とても低い温度では粗くなることが、教科書で知られています（本記事による解説。Allen と Tildesley の教科書の一般的な説明にもとづく）［参考資料6］。電子の状態を毎回計算する第一原理MD（後述）という方法もありますが、計算は重くなります。
3. 力はどこから来るか：エネルギーの坂
では、「いま各原子にかかる力」は、どうやって決めるのでしょうか。考え方は、坂の上の玉と同じです。玉は、坂を下る向きに力を受けます。原子も、エネルギー（ポテンシャルエネルギー）の坂を下る向きに力を受けます。つまり、力は「エネルギーの坂の傾き」に、向きを逆にしたものです。原子の位置を変えたときのエネルギーを表す関数をポテンシャルと呼びます。
いちばん単純な例が、2つの原子のあいだの力です。2つの原子は、近すぎると強く押し合い、離れていると弱く引き合います。その中間に、エネルギーがいちばん低くなる「つり合い」の距離があります。この形を簡単な式で表したものが、レナード＝ジョーンズ型のポテンシャルです（式の原型は Jones、1924）事実［参考資料12］。
- 近づくと：r＝σ のとき U＝0 ですが、さらに近づくとエネルギーは急に高くなります。力は、離す向きに大きくなります（r＝σ で 24ε/σ）本稿の計算
- つり合いの距離：r≒1.12σ（＝2の6乗根×σ）で、エネルギーはいちばん低い −ε になり、力はゼロになります本稿の計算
- 離れると：引力は弱くなります。r＝1.5σ で力は −1.16（引き合う向き）、r＝2.5σ では U は −0.016ε で、底の深さの約1.6%です本稿の計算。MDでは、この遠くの原子との力を省略（打ち切り）して計算を速くすることがよくあります。本稿の簡易計算も 2.5σ で打ち切りました
σ は「原子の大きさの目安になる距離」、ε は「エネルギーの底の深さ（結びつきの強さの目安）」です。物質ごとに違う値を当てはめて使います。図では、どちらも「1」とした相対値で描いています。本物の金属や高分子の原子は、1対ずつの力の足し合わせでは表せない部分があり（次の章）、このレナード＝ジョーンズ型は考え方を学ぶための最も単純な例と思ってください（本記事による解説）。
4. 時間刻み：なぜ約1フェムト秒なのか
図2の「ごく短い時間 Δt」は、どのくらい短くすればよいのでしょうか。大きくすれば、同じ計算で長い時間を動かせます。しかし、大きくしすぎると、計算がおかしくなるのです。
理由は、原子がとても速く振動しているからです。たとえば、C–H 結合の伸び縮みの波数は、教科書的な目安で約 3000 cm⁻¹ です。周期は 1/(光速×波数) で、約 11 fs になります本稿の計算。1周期を10コマ程度に分けないと、振動の動きを追えません。1 fs 前後という刻みは、こうした速い振動に合わせたものです。LAMMPS の既定値も、real 単位で 1.0 fs、metal 単位で 0.001 ps です事実［参考資料9］。
- 刻みを細かくすると：Δt＝T₀/50（T₀ は振動の周期）なら、200周期にわたって、エネルギーのずれは最大0.22%でした本稿の計算。刻みを半分にすると、ずれは約4分の1になります（0.55×10⁻³ → 0.14×10⁻³）本稿の計算
- 刻みを粗くすると：Δt＝T₀/8 では、ずれは最大13%になり、Δt＝T₀/6 以上では計算が破綻しました本稿の計算
- 方法による違い：速度ベルレ法は、ずれが大きくなっても、一定の幅で振れるだけで、時間とともに増え続けません。これに対し、最も素朴な更新法（単純なオイラー法）は、同じ Δt＝T₀/50 でも数周期でエネルギーが急増しました本稿の計算
外から熱や仕事を加えない孤立した系では、全エネルギー（運動エネルギー＋ポテンシャルエネルギー）は、一定に保たれるはずです。計算でこれが保たれなくなっていれば、数値計算の誤差が大きすぎる目印になります。MDの結果をグラフにする前に、エネルギーが保存されているかを確かめるのは、基本の点検です（本記事による解説）。
水素を含む結合の長さを固定する方法（SHAKE など）を使うと、速い振動を計算から外せるため、刻みを長くできる場合があります事実［参考資料18］。ただし、固定した分、その振動は表せなくなります（本記事による解説）。
5. 力場の役割と種類
MDの結果の良し悪しは、力の決め方で決まります。力の決め方（原子のあいだのエネルギーを表す式と、その定数のセット）を力場（りきば、フォースフィールド）と呼びます。同じ MD のソフトウェアでも、力場が違えば、まったく違う結果が出ます。
欧州の規格化団体 CEN が2018年に公開した文書（CEN Workshop Agreement 17284）は、この関係を次のように整理しています。運動方程式は「物理の式」で、どんな材料にも共通です。一方、力場は「材料に固有の関係式」で、その式を作った材料の種類に対してだけ有効です事実［参考資料34］。力場は、原子の組み合わせ・温度・圧力・構造など、作られた条件の外では信頼できないと考えるのが基本です（本記事による解説）。なお、この文書は、正式な規格ではなく、ワークショップの合意文書です［参考資料34］。
| 方法 | 考え方 | 代表的な使い先（原論文の題名・内容による） | 読み方の注意 |
|---|---|---|---|
| レナード＝ジョーンズ型 | 2原子のあいだに、近いと反発・離れると引力。式は簡単 | 気体の状態方程式の解析（Jones、1924）。希ガスの液体のモデル（Rahman、1964）［参考資料12・2］事実 | 考え方を学ぶ・基本の性質を調べる用。実際の金属や共有結合の物質をそのまま表すものではない（本記事による解説） |
| EAM（埋め込み原子法） | 原子が、周りの電子の密度の中に「埋め込まれている」と考える | 金属の不純物・表面・欠陥の解析（Daw と Baskes、1984）［参考資料13］事実 | 金属向け。対象の元素・合金の組成ごとに作られている |
| スティリンジャー＝ウェーバー、Tersoff | 原子の周りの結合の様子（角度など）も式に入れる | シリコンの局所的な秩序（Stillinger と Weber、1985）、共有結合の物質の構造とエネルギー（Tersoff、1988）［参考資料14・15］事実 | 共有結合の結晶向け。適用できる元素・構造は、原論文で確認する |
| ReaxFF（反応力場） | 結合の強さが距離によって変わる式で、結合の切れ・生成を表す | 炭化水素の反応（van Duin ら、2001）［参考資料16］事実 | 反応を扱えるが、式と定数が複雑で、対象の系ごとに検証が必要（本記事による解説） |
| 第一原理MD | 各ステップで電子の状態を計算して力を出す。Car と Parrinello（1985）が、MDと密度汎関数理論を組み合わせた方法を提案 | 電子の働きが本質の系。別記事（密度汎関数理論の基礎）［参考資料17］事実 | 最も重く、扱える原子数と時間は小さい（本記事による解説） |
| 機械学習ポテンシャル | 第一原理計算のデータから、力とエネルギーの関係を学習 | 別記事（機械学習原子間ポテンシャル） | 学習データの範囲の外は不確か。別記事で詳しく扱う |
各方法の式の詳細は、原論文を確認してください。本稿は、原論文の題名と内容の範囲で整理しています。
力場を選ぶときに頼れる公開の場もあります。NIST の「原子間ポテンシャルのリポジトリ」は、ポテンシャル（力場）と関連ファイル、その質と適用範囲を判断するための評価ツールを公開しています事実［参考資料30］。OpenKIM は、機械学習によるものを含む原子間ポテンシャルを集めた、審査つきのリポジトリです事実［参考資料31］。力場の選び方を述べた論文もあります［参考資料29］。
どの力場も、作られた対象（元素の組み合わせ、相、温度の範囲）で、何を再現するように調整されたかが決まっています。CEN の文書が述べるとおり、材料に固有の関係式は、その材料の種類に対してだけ有効です［参考資料34］。ある合金の融点と弾性率を再現する力場でも、界面や、不純物を多く含む試料で正しいとは限りません（本記事による解説）。選ぶときは、自分の材料と条件で検証されているかを、原論文や NIST・OpenKIM の評価で確かめます。
6. 温度・圧力を保つ方法と、小さな箱の工夫
現実の実験は、「室温で」「大気圧で」行われます。MDでも、箱の中の温度や圧力を一定に保って計算したいことが多くあります。何を一定にするかの組合せをアンサンブルと呼びます。粒子数 N・体積 V・エネルギー E を一定にするものが NVE（孤立系）、N・V・温度 T を一定にするものが NVT、N・圧力 P・T を一定にするものが NPT です（本記事による解説）［参考資料6］。
- 温度を一定にする：系を「熱浴」とつなぎ、原子の速度を調整するサーモスタットを使います。代表的な方法に、Nosé（1984）と Hoover（1985）の方法があります事実［参考資料19・20］
- 圧力を一定にする：箱の大きさを少しずつ変えるバロスタットを使います（本記事による解説）［参考資料6］
- 小さな箱をくり返す：計算できる原子は数千〜数百万個で、現実の材料よりはるかに小さな塊です。そこで、箱の反対側の面をつなげ、同じ箱が無限にくり返されている、と見なす周期境界条件を使い、表面の影響を消します。ただし、箱が小さいと、結果が箱の大きさに左右されます（本記事による解説）［参考資料6］
7. 時間と大きさの限界
MDの最大の弱点は、届く時間と大きさが限られることです。1 fs のステップで 1 ns（ナノ秒）を動かすには 10⁶ 回、1 μs（マイクロ秒）なら 10⁹ 回、1 ms なら 10¹² 回、1 秒なら 10¹⁵ 回の計算が要ります本稿の計算。大きさも同じです。鉄 1 mm³ には、約 8.5×10¹⁹ 個の原子があります本稿の計算。
- 時間の記録：専用に作られた計算機を使うと、はるかに長い時間に届きます。Shaw らの論文（2010）は、専用機で、タンパク質 BPTI の1 ミリ秒の全原子MD（古典力場）を報告しています事実［参考資料21］。1 fs のステップなら約 10¹² 回の計算に相当します（本稿の換算）本稿の計算
- 大きさの記録：Germann と Kadau（2008）は、BlueGene/L（212,992 個のプロセッサ）で、1兆個の原子（レナード＝ジョーンズ粒子）のMDを行ったと報告しています事実［参考資料22］。しかし、1兆個の原子は、鉄なら一辺が約2.3 μm の立方体にすぎません本稿の計算
- 実用で効く時間との差：拡散、クリープ、腐食、疲労のように、秒〜年で進む現象は、ふつうのMDでは直接届きません（本記事による解説）。ごくまれにしか起きない出来事（原子の跳躍、核生成）は、シミュレーション中に一度も起きないことがあります
この壁を越える工夫として、まれな出来事を加速する方法があります。たとえば、ハイパーダイナミクス（Voter、1997）は、まれにしか起きない出来事の進行を加速する方法で、メタダイナミクス（Laio と Parrinello、2002）は、エネルギーの谷を「埋めて」脱出しやすくし、自由エネルギーの形を調べる方法です事実［参考資料23・24］。ただし、どの変数で加速するかなど、使う側の選択が結果に影響します（本記事による解説）。
計算機が速くなっても、ステップ数が10⁶倍になれば、計算も10⁶倍かかります（本稿の計算の考え方）。1 ms に達した例は、専用の計算機を使った特別な記録です［参考資料21］。固体の中の遅い拡散や腐食のように、秒以上かかる現象をそのまま追うのは、通常の方法では難しい、と考えるのが現実的です。高い温度で計算して速める、加速法を使う、別の方法（マルチスケール、別記事）と組み合わせる、といった工夫と、その近似の検証が要ります（本記事による解説）。
8. MDから何が分かるか
MDで得られるのは、全原子の位置と速度の時間変化です。そこから、次のような量を計算できます。ここでは、本稿で行った簡単な計算（256個の粒子）で、液体と固体の違いを見てみます。
- 原子の並び方（動径分布関数）：ある原子から距離 r のところに、別の原子がどれだけいるかを表します。固体では山が高く、遠くまで規則的に並びます（最大の山 約5.9）。液体では、近くの山だけが残り、遠くでは1に近づきます（最大の山 約3.0）本稿の計算
- 動きやすさ（平均二乗変位と拡散係数）：原子が動いた距離の2乗の平均です。液体では、時間に比例して増えます。3次元では「拡散係数＝傾き÷6」と置く関係が、教科書に載っています［参考資料6］。本稿の簡易計算では、液体で約0.026（規格化した単位）でした。固体では、ほとんど増えません（5000ステップで約0.011）本稿の計算。短時間・小さな系の目安の値で、誤差の評価はしていません
| 調べたい量 | MDでの見方 | 原論文・総説での例 |
|---|---|---|
| 拡散・イオンの動き | MSD の傾きから拡散係数。軌跡から、イオンが動く道筋を調べる | 固体電解質で、イオンが複数個で協調して動くことを、第一原理に基づく計算で示した研究（He ら、2017）［参考資料27］事実 |
| 熱伝導 | 熱の流れを計算する方法が複数あり、結果が方法に左右される | 原子レベルの計算で熱伝導率を求める方法の比較（Schelling ら、2002）［参考資料25］事実 |
| 変形・強度 | 箱を引き伸ばし、応力とひずみの関係、転位などの欠陥の動きを見る | 転位のからみ合い（多重接合）の存在を、転位動力学と原子シミュレーションで予測し、モリブデンの電子顕微鏡観察で確認（Bulatov ら、2006）［参考資料26］事実 |
| 高分子の動き | 鎖を粗い玉とバネで表し（粗視化）、絡み合いを調べる | 絡み合った高分子溶融体の動き（Kremer と Grest、1990）［参考資料28］事実 |
| 生体分子の動き | タンパク質の構造の揺らぎ・折りたたみ | McCammon ら（1977）、Karplus と McCammon の総説（2002）［参考資料3・4］事実 |
Bulatov らの例は、計算で予測した現象を、実験で確かめた例です。多くのMDの結果は、実験との直接の比較がなく、計算上の予測にとどまります（本記事による解説）。
固体は、原子が決まった位置に並ぶので、「ある距離に原子がいる」確率に、はっきりした山と谷ができます。液体は、近くの原子とは一定の距離を保ちますが、少し離れると、位置がばらけて山が消えます（図7左）。X線回折の測定でも、同じ種類の情報（原子の並び方）が得られるので、計算と実測を比べやすい量の一つです（本記事による解説）。
9. 研究の進め方と、確かめるところ
MDの研究は、計算を回すだけではなく、確認が多い仕事です。LAMMPS（Thompson ら、2022）や GROMACS（Abraham ら、2015）は、広く使われる計算ソフトウェアとして、論文で公開されています事実［参考資料10・11］。ソフトウェアを使えば計算はできますが、結果を信じてよいかは、使う側の確認にかかっています。
- 初心者向けの指針：分子シミュレーションを始める人が、始める前と最中に知っておくべき点をまとめた文書と、チェックリストが公開されています（Braun ら）事実［参考資料32］
- 統計の確認：不確かさと、サンプリングの質（動かした時間が十分か）の評価は、分子シミュレーションで不可欠だとする指針があります（Grossfield ら）。この指針は、力場の不正確さのような系統的な誤差は扱っていないと明記しています事実［参考資料33］。つまり、誤差棒が小さくても、力場が合っているとは限りません
- 計算の記録：CEN の文書は、シミュレーションの記録の様式（MODA）を示しています。使い道（ユースケース）、モデル、時間刻み、ソフトウェアなどを項目として書き残す形です事実［参考資料34］
10. 素材技術者から見て
MDで直接出るのは、原子の位置と速度です。ここから、拡散係数、密度、弾性率、熱伝導率、ガラス転移温度の目安などを計算で導きます。自社の実験で測っている量と、MDで計算できる量が一対一で対応するかを、計算を始める前に表にします。対応しない量（たとえば、観測できない原子の動き）は、「仮説」として扱い、実験で確かめる設計（どの試料・条件で測れば検証できるか）を考えておきます（本記事による解説）。
公開されている力場は、純金属や単純な合金、標準的な高分子のような代表的な系で作られていることが多いと考えられます（推測。個別の確認はしていません）。自社の材料には、添加物、不純物、界面、複数の相が入ります。力場を選ぶときは、どの元素の組合せ・温度・構造で、何を再現するように作られたかを原論文で確認し、NIST や OpenKIM の評価を見ます［参考資料30・31］。そのうえで、自社の簡単な実測値（密度、弾性率、融点など）を再現するかを、まず試します。再現しない力場の拡散係数は、信頼しません。
実験の時間（秒〜年）に届かないとき、「温度を上げて速く進める」という工夫が考えられます。ただし、高温では、別の機構が働くことがあります（本記事による解説）。温度を変えて、何本かの計算をして傾向を見る、加速法と組み合わせる、といった設計をし、結果を実験の温度へ外挿するときに、機構が変わっていないかを、実測との比較で確認します。
MDの結果は、初期の構造（結晶か、ガラス状か、結晶粒の大きさ、不純物の位置）に影響されます。自社の試料について、製法、熱履歴、組成の分析値、観察像を、実験ノートに記録しておくと、計算の初期構造を、現実に近づけられます。うまくいかなかった計算や実験の条件も、失敗データとして残してください。力場の妥当性や、適用範囲の判断に役立ちます。
11. 限界と、確認できなかったこと／まとめ
- MDは、原子をニュートンの運動方程式で、約1 fs ずつ動かす計算です。速度ベルレ法のような積分法を使い、刻みが粗すぎるとエネルギーが保存されず、計算が破綻します事実本稿の計算
- 結果は、力場で決まります。力場は、作られた材料の種類に対してだけ有効です。自社の材料と条件で検証されているかを、確かめます事実
- 時間と大きさに限界があります。1 ms は 10¹² ステップ、1兆個の原子は鉄で一辺約2.3 μm です。実験の時間（秒〜年）には、通常は届きません本稿の計算
- MDで分かる量（拡散、熱伝導、変形、構造）は、実測と比べて初めて意味が出ます。計算上の予測と、実験で確かめた結果を区別します事実
確認できなかったこと
- Allen と Tildesley の教科書の個々の記述（古典近似が粗くなる条件、アンサンブル、周期境界条件、拡散係数の関係式）は、書籍の本文を本稿の作業では確認していません。一般的な説明として書いており、必要なときは原典で確認してください
- C–H 結合の波数（約3000 cm⁻¹）は、教科書的な目安で、個別の文献は確認していません。周期の計算は、この値からの換算です
- 簡易MDの結果（g(r) の山の高さ、拡散係数など）は、短時間・小さな系の目安で、文献の値との比較、誤差の評価、箱の大きさ依存の確認はしていません
- 実際の素材（自社の合金・高分子）で、どの力場がどの程度の精度かは、本稿では一般論としても確認していません
- MDの力場や計算の手順そのものを定める、国際規格は確認できませんでした
12. 評価・標準の読み方
MDの「ものさし」は、まず計算が数値として破綻していないか（エネルギーの保存、平衡に達したか）、次に統計の信頼性（動かした時間は十分か、複数回やっても同じか）、そして力場の妥当性（既知の実測を再現するか）です。最後に、実測と同じ量で比べられるかです。確認日は2026年10月22日です。
| 評価項目 | 代表的なベンチマーク・評価手法（名称・出典） | 標準・指針（ISO/IEC・NIST・CEN など） | 読み方の注意 |
|---|---|---|---|
| 数値計算の健全さ（エネルギー保存・時間刻み） | 孤立系のエネルギーが保たれるかの確認（本稿の図4）。LAMMPS の既定の刻み：real 単位で 1.0 fs、metal 単位で 0.001 ps［参考資料9］本稿の計算事実 | 確認できなかった（時間刻みを定める規格は確認していない） | 既定値は出発点で、系ごとに確かめる。軽い原子を含む系は短くする必要がある場合がある |
| 力場・ポテンシャルの質 | NIST 原子間ポテンシャルのリポジトリ：評価ツールつき［参考資料30］。OpenKIM：審査つきの公開リポジトリ［参考資料31］事実 | 確認できなかった（力場の認証制度は確認していない）。CEN CWA 17284 は、力場を「材料に固有の関係式」と位置づける用語の整理［参考資料34］ | 公開されていても、自分の材料・条件で検証されたとは限らない |
| 統計の信頼性（サンプリング・不確かさ） | Grossfield ら：不確かさとサンプリング品質の評価の指針［参考資料33］事実 | 確認できなかった（国際規格は確認していない） | 誤差棒は、力場の系統的な誤差を含まない（同指針が明記） |
| 計算の準備・再現性 | Braun ら：分子シミュレーションの基礎の指針とチェックリスト［参考資料32］事実。CEN CWA 17284：シミュレーションの記録様式 MODA［参考資料34］事実 | CWA 17284:2018（CEN Workshop Agreement、2018年4月）。正式な規格ではないとされる［参考資料34］事実 | MODA は記録の様式で、結果の正しさを保証するものではない |
| 実測との比較 | 同じ量（密度、弾性率、拡散係数、X線回折など）を実測と比べる（解説） | 確認できなかった | 観測できない量は、仮説として扱う。計算と実験の差の原因（力場、系の大きさ、時間）を切り分ける |
表の整理は本稿の計算（整理）。番号・題名・状態は公式ページ・原論文で確認（確認日 2026年10月22日）。CWA の本文は公開されており、題名・適用範囲・用語の定義を確認しました。ASTM の規格は、本稿の範囲では該当を確認していません。
読み方の注意（素材技術者向け）
- 「誤差棒が小さい」と「合っている」は別です：統計の誤差は、動かす時間を伸ばせば小さくなります。しかし、力場の誤りによる系統的な誤差は残ります［参考資料33］
- 公開された力場のランキングは、自社の材料では当てになりません：各力場は、作られた対象の性質を再現するように調整されています。自社の材料の簡単な実測値で、まず試します
- 計算上の予測と、実験で確かめた結果を区別してください：MDの多くの結果は、実験との直接の比較がない計算上の予測です。実験で確かめた例（Bulatov ら）は、そのことが明記されます［参考資料26］
- 標準や指針は、推奨か要求か：本稿で挙げた指針・CWAは、推奨・用語整理・記録様式で、法的な拘束力をもつ規格ではありません［参考資料34］
13. 用語集
- 分子動力学（MD）
- 原子を、ニュートンの運動方程式で時間を追って動かす計算。molecular dynamics。
- ニュートンの運動方程式
- 力＝質量×加速度。力が決まれば、加速度が決まる。
- ポテンシャル（ポテンシャルエネルギー）
- 原子の位置から決まるエネルギー。その坂（傾き）から力が決まる。
- 力場（フォースフィールド）
- 原子のあいだのエネルギーを表す式と、その定数のセット。作られた材料の種類に対してだけ有効。
- レナード＝ジョーンズ型
- 近いと反発し、離れると弱く引き合う、最も単純なポテンシャルの形。
- EAM（埋め込み原子法）
- 原子が周りの電子密度の中に埋め込まれていると考える、金属向けの方法。
- ReaxFF
- 結合の強さが距離で変わる式で、結合の切れ・生成を表す反応力場。
- 第一原理MD
- 各ステップで電子の状態を計算して力を出すMD。重いが、電子の働きに近い。
- 時間刻み（Δt）
- 1ステップで進める時間。約1 fs。粗すぎると計算が破綻する。
- フェムト秒（fs）
- 10⁻¹⁵秒。1 ps（ピコ秒）＝1000 fs、1 ns（ナノ秒）＝10⁶ fs。
- 速度ベルレ法
- 運動方程式を少しずつ進める標準的な方法。エネルギーのずれが増え続けにくい。
- アンサンブル（NVE・NVT・NPT）
- 何を一定にして計算するかの組合せ（粒子数・体積・エネルギー・温度・圧力）。
- サーモスタット／バロスタット
- 温度／圧力を一定に保つ工夫。
- 周期境界条件
- 箱がくり返されていると見なし、表面の影響を消す工夫。
- 打ち切り（カットオフ）
- 遠くの原子との力を省略して、計算を速くする工夫。
- 動径分布関数 g(r)
- ある原子から距離 r に別の原子がいる確率（相対値）。原子の並び方を表す。
- 平均二乗変位（MSD）
- 原子が動いた距離の2乗の平均。拡散のしやすさの指標。
- 拡散係数
- 原子の動きやすさ。3次元では MSD の傾き÷6。
- 稀事象
- めったに起きない出来事（原子の跳躍、核生成など）。直接のMDでは捉えにくい。
- メタダイナミクス／ハイパーダイナミクス
- まれな出来事を加速するMDの工夫。
- 粗視化
- 複数の原子を1つの玉にまとめて、大きな系を扱う方法。
- MODA
- materials modelling data。CENの文書が示す、シミュレーションの記録様式。
14. 参考資料（一次情報）
- Alder B. J., Wainwright T. E. "Phase Transition for a Hard Sphere System", J. Chem. Phys. 27, 1208–1209 (1957), DOI 10.1063/1.1743957 — doi.org/10.1063/1.1743957
- Rahman A. "Correlations in the Motion of Atoms in Liquid Argon", Phys. Rev. 136, A405–A411 (1964), DOI 10.1103/PhysRev.136.A405 — doi.org/10.1103/PhysRev.136.A405
- McCammon J. A., Gelin B. R., Karplus M. "Dynamics of folded proteins", Nature 267, 585–590 (1977), DOI 10.1038/267585a0 — doi.org/10.1038/267585a0
- Karplus M., McCammon J. A. "Molecular dynamics simulations of biomolecules", Nat. Struct. Biol. 9, 646–652 (2002), DOI 10.1038/nsb0902-646 — doi.org/10.1038/nsb0902-646
- The Royal Swedish Academy of Sciences "The Nobel Prize in Chemistry 2013, press release (9 October 2013)" — www.nobelprize.org/prizes/chemistry/2013/press-release/
- Allen M. P., Tildesley D. J. "Computer Simulation of Liquids, 2nd edition", Oxford University Press (2017), DOI 10.1093/oso/9780198803195.001.0001 — doi.org/10.1093/oso/9780198803195.001.0001
- Verlet L. "Computer "Experiments" on Classical Fluids. I. Thermodynamical Properties of Lennard-Jones Molecules", Phys. Rev. 159, 98–103 (1967), DOI 10.1103/PhysRev.159.98 — doi.org/10.1103/PhysRev.159.98
- LAMMPS documentation "run_style command (verlet = velocity form of the Stoermer-Verlet algorithm)" — docs.lammps.org/run_style.html
- LAMMPS documentation "timestep command (default 0.001 ps for metal units, 1.0 fs for real units)" — docs.lammps.org/timestep.html
- Thompson A. P. ほか "LAMMPS - a flexible simulation tool for particle-based materials modeling at the atomic, meso, and continuum scales", Comput. Phys. Commun. 271, 108171 (2022), DOI 10.1016/j.cpc.2021.108171 — doi.org/10.1016/j.cpc.2021.108171
- Abraham M. J. ほか "GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers", SoftwareX 1–2, 19–25 (2015), DOI 10.1016/j.softx.2015.06.001 — doi.org/10.1016/j.softx.2015.06.001
- Jones J. E. "On the determination of molecular fields. II. From the equation of state of a gas", Proc. R. Soc. Lond. A 106, 463–477 (1924), DOI 10.1098/rspa.1924.0082 — doi.org/10.1098/rspa.1924.0082
- Daw M. S., Baskes M. I. "Embedded-atom method: Derivation and application to impurities, surfaces, and other defects in metals", Phys. Rev. B 29, 6443–6453 (1984), DOI 10.1103/PhysRevB.29.6443 — doi.org/10.1103/PhysRevB.29.6443
- Stillinger F. H., Weber T. A. "Computer simulation of local order in condensed phases of silicon", Phys. Rev. B 31, 5262–5271 (1985), DOI 10.1103/PhysRevB.31.5262 — doi.org/10.1103/PhysRevB.31.5262
- Tersoff J. "New empirical approach for the structure and energy of covalent systems", Phys. Rev. B 37, 6991–7000 (1988), DOI 10.1103/PhysRevB.37.6991 — doi.org/10.1103/PhysRevB.37.6991
- van Duin A. C. T. ほか "ReaxFF: A Reactive Force Field for Hydrocarbons", J. Phys. Chem. A 105, 9396–9409 (2001), DOI 10.1021/jp004368u — doi.org/10.1021/jp004368u
- Car R., Parrinello M. "Unified Approach for Molecular Dynamics and Density-Functional Theory", Phys. Rev. Lett. 55, 2471–2474 (1985), DOI 10.1103/PhysRevLett.55.2471 — doi.org/10.1103/PhysRevLett.55.2471
- Ryckaert J.-P., Ciccotti G., Berendsen H. J. C. "Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes", J. Comput. Phys. 23, 327–341 (1977), DOI 10.1016/0021-9991(77)90098-5 — doi.org/10.1016/0021-9991(77)90098-5
- Nosé S. "A unified formulation of the constant temperature molecular dynamics methods", J. Chem. Phys. 81, 511–519 (1984), DOI 10.1063/1.447334 — doi.org/10.1063/1.447334
- Hoover W. G. "Canonical dynamics: Equilibrium phase-space distributions", Phys. Rev. A 31, 1695–1697 (1985), DOI 10.1103/PhysRevA.31.1695 — doi.org/10.1103/PhysRevA.31.1695
- Shaw D. E. ほか "Atomic-Level Characterization of the Structural Dynamics of Proteins", Science 330, 341–346 (2010), DOI 10.1126/science.1187409 — doi.org/10.1126/science.1187409
- Germann T. C., Kadau K. "Trillion-atom molecular dynamics becomes a reality", Int. J. Mod. Phys. C 19, 1315–1319 (2008), DOI 10.1142/S0129183108012911 — doi.org/10.1142/S0129183108012911
- Voter A. F. "Hyperdynamics: Accelerated Molecular Dynamics of Infrequent Events", Phys. Rev. Lett. 78, 3908–3911 (1997), DOI 10.1103/PhysRevLett.78.3908 — doi.org/10.1103/PhysRevLett.78.3908
- Laio A., Parrinello M. "Escaping free-energy minima", Proc. Natl. Acad. Sci. USA 99, 12562–12566 (2002), DOI 10.1073/pnas.202427399 — doi.org/10.1073/pnas.202427399
- Schelling P. K., Phillpot S. R., Keblinski P. "Comparison of atomic-level simulation methods for computing thermal conductivity", Phys. Rev. B 65, 144306 (2002), DOI 10.1103/PhysRevB.65.144306 — doi.org/10.1103/PhysRevB.65.144306
- Bulatov V. V. ほか "Dislocation multi-junctions and strain hardening", Nature 440, 1174–1178 (2006), DOI 10.1038/nature04658 — doi.org/10.1038/nature04658
- He X., Zhu Y., Mo Y. "Origin of fast ion diffusion in super-ionic conductors", Nat. Commun. 8, 15893 (2017), DOI 10.1038/ncomms15893 — doi.org/10.1038/ncomms15893
- Kremer K., Grest G. S. "Dynamics of entangled linear polymer melts: A molecular-dynamics simulation", J. Chem. Phys. 92, 5057–5086 (1990), DOI 10.1063/1.458541 — doi.org/10.1063/1.458541
- Becker C. A., Tavazza F., Trautt Z. T., Buarque de Macedo R. A. "Considerations for choosing and using force fields and interatomic potentials in materials science and engineering", Curr. Opin. Solid State Mater. Sci. 17, 277–283 (2013), DOI 10.1016/j.cossms.2013.10.001 — doi.org/10.1016/j.cossms.2013.10.001
- NIST (National Institute of Standards and Technology) "Interatomic Potentials Repository" — www.ctcms.nist.gov/potentials/
- OpenKIM "Knowledgebase of Interatomic Models (OpenKIM)" — openkim.org/
- Braun E. ほか "Best Practices for Foundations in Molecular Simulations [Article v1.0]", Living J. Comput. Mol. Sci. 1, 5957, DOI 10.33011/livecoms.1.1.5957 — doi.org/10.33011/livecoms.1.1.5957
- Grossfield A. ほか "Best Practices for Quantification of Uncertainty and Sampling Quality in Molecular Simulations [Article v1.0]", Living J. Comput. Mol. Sci. 1, 5067, DOI 10.33011/livecoms.1.1.5067 — doi.org/10.33011/livecoms.1.1.5067
- CEN (European Committee for Standardization) "CEN Workshop Agreement CWA 17284:2018, Materials modelling - Terminology, classification and metadata" (2018) — www.cencenelec.eu/media/CEN-CENELEC/CWAs/RI/cwa17284_2018.pdf
15. 主張と出典の対応表（監査）
| 本文の主張 | 根拠 | 区分 |
|---|---|---|
| MD の1ステップの既定値：LAMMPS の timestep の既定値は、metal 単位で 0.001 ps、real 単位で 1.0 fs | 出典9 https://docs.lammps.org/timestep.html | 事実 |
| LAMMPS の既定の積分法 verlet は、Stoermer-Verlet 法の速度版（速度ベルレ法） | 出典8 https://docs.lammps.org/run_style.html | 事実 |
| ベルレ法は Verlet（1967）の論文に由来する。初期の MD の例：Alder と Wainwright（1957）の硬い球の系、Rahman（1964）の液体アルゴン、McCammon ら（1977）のタンパク質の動力学 | 出典7 https://doi.org/10.1103/PhysRev.159.98 出典1 https://doi.org/10.1063/1.1743957 出典2 https://doi.org/10.1103/PhysRev.136.A405 出典3 https://doi.org/10.1038/267585a0 | 事実 |
| McCammon らの研究は、経験的なポテンシャルエネルギー関数で、折りたたまれたタンパク質（BPTI）の運動方程式を解いた。Karplus と McCammon の総説：生体分子の MD の歴史と使われ方 | 出典3 https://doi.org/10.1038/267585a0 出典4 https://doi.org/10.1038/nsb0902-646 | 事実 |
| 2013年のノーベル化学賞：「複雑な化学系のための多重スケールモデルの開発」。古典物理と量子物理を並べて使えるようにした（本稿では、MD の力を電子の状態に近づける方向の説明で言及） | 出典5 https://www.nobelprize.org/prizes/chemistry/2013/press-release/ | 事実 |
| レナード＝ジョーンズ型のポテンシャルの式の原型：Jones（1924）、気体の状態方程式からの分子間力の決定 | 出典12 https://doi.org/10.1098/rspa.1924.0082 | 事実 |
| 図3の数値：U＝0（r＝σ）、F＝24ε/σ（r＝σ）、最小 −ε（r≒1.12σ）、r＝1.5σ で F＝−1.16、r＝2.5σ で U＝−0.016ε。式 U＝4ε[(σ/r)¹²−(σ/r)⁶] から本稿で計算 | 本稿の計算 | 本稿の計算 |
| C–H 伸縮の周期約11 fs：波数 約3000 cm⁻¹（教科書的な目安）から、周期＝1/(光速×波数) で本稿が換算。個別の文献は確認していない | 本稿の計算（教科書的な目安からの換算） | 本稿の計算 |
| 図4の数値（2原子のレナード＝ジョーンズ型、初期1.3σ、200周期）：Δt＝T₀/50 で最大0.22%、T₀/100 で0.055%、T₀/200 で0.014%、T₀/8 で13%、T₀/6 以上で破綻。単純なオイラー法は数周期でエネルギーが急増。速度ベルレ法は一定の幅で振れる | 本稿の計算（速度ベルレ法の数値積分） | 本稿の計算 |
| 水素を含む結合の長さを固定する方法（SHAKE など）：Ryckaert らの拘束つき運動方程式の積分法。刻みを長くできる場合があるという記述は解説 | 出典18 https://doi.org/10.1016/0021-9991(77)90098-5 | 事実 |
| CEN Workshop Agreement 17284:2018（2018年4月）：運動方程式は「物理の式」、力場は「材料に固有の関係式」で、その材料の種類に対してだけ有効。正式な規格ではない。MODA はシミュレーションの記録様式で、使い道・モデル・ソルバー・後処理を記す。Verlet 積分は時間の離散化のみを要するソルバーの例 | 出典34 https://www.cencenelec.eu/media/CEN-CENELEC/CWAs/RI/cwa17284_2018.pdf | 事実 |
| 各力場の原論文：EAM（Daw と Baskes、1984）、シリコン（Stillinger と Weber、1985）、共有結合（Tersoff、1988）、ReaxFF（van Duin ら、2001）、第一原理MD（Car と Parrinello、1985） | 出典13 https://doi.org/10.1103/PhysRevB.29.6443 出典14 https://doi.org/10.1103/PhysRevB.31.5262 出典15 https://doi.org/10.1103/PhysRevB.37.6991 出典16 https://doi.org/10.1021/jp004368u 出典17 https://doi.org/10.1103/PhysRevLett.55.2471 | 事実 |
| NIST 原子間ポテンシャルのリポジトリは、ポテンシャルと関連ファイル、質と適用範囲を判断する評価ツールを公開。OpenKIM は機械学習を含む原子間ポテンシャルの審査つきリポジトリ。力場の選び方を述べた Becker らの論文 | 出典30 https://www.ctcms.nist.gov/potentials/ 出典31 https://openkim.org/ 出典29 https://doi.org/10.1016/j.cossms.2013.10.001 | 事実 |
| サーモスタット：Nosé（1984）、Hoover（1985）。アンサンブル・バロスタット・周期境界条件の一般的な説明は、Allen と Tildesley の教科書にもとづく解説（書籍の本文は確認していない） | 出典19 https://doi.org/10.1063/1.447334 出典20 https://doi.org/10.1103/PhysRevA.31.1695 出典6 https://doi.org/10.1093/oso/9780198803195.001.0001 | 事実 |
| 1 ns＝10⁶ステップ、1 μs＝10⁹、1 ms＝10¹²、1 s＝10¹⁵（1 fs のステップ）。鉄 1 mm³ は約8.5×10¹⁹ 個、10¹² 個は一辺約2.3 μm の立方体（密度 7.87 g/cm³、モル質量 55.85 g/mol からの概算） | 本稿の計算 | 本稿の計算 |
| Shaw ら：専用機での全原子MD。BPTI の 1 ミリ秒の MD。WW ドメインの折りたたみ・展開の観察 | 出典21 https://doi.org/10.1126/science.1187409 | 事実 |
| Germann と Kadau：BlueGene/L（212,992 個の PowerPC440 プロセッサ）、SPaSM、1兆個の原子のレナード＝ジョーンズ粒子のMD | 出典22 https://doi.org/10.1142/S0129183108012911 | 事実 |
| まれな出来事の加速：ハイパーダイナミクス（Voter、1997）、メタダイナミクス（Laio と Parrinello、2002） | 出典23 https://doi.org/10.1103/PhysRevLett.78.3908 出典24 https://doi.org/10.1073/pnas.202427399 | 事実 |
| 図7の簡易MD（256粒子、周期境界、打ち切り2.5σ、Δt＝0.004）：g(r) の最大の山は固体で約5.9、液体で約3.0。液体の拡散係数は約0.026、固体の MSD は5000ステップで約0.011。短時間・小さな系の目安で、誤差評価なし | 本稿の計算（簡易MD） | 本稿の計算 |
| 固体電解質でのイオンの協調的な動き（He ら、2017）。熱伝導率を求める方法の比較（Schelling ら、2002）。転位多重接合：転位動力学と原子シミュレーションで予測し、モリブデンの TEM で確認（Bulatov ら、2006）。絡み合った高分子溶融体（Kremer と Grest、1990） | 出典27 https://doi.org/10.1038/ncomms15893 出典25 https://doi.org/10.1103/PhysRevB.65.144306 出典26 https://doi.org/10.1038/nature04658 出典28 https://doi.org/10.1063/1.458541 | 事実 |
| LAMMPS（Thompson ら、2022）、GROMACS（Abraham ら、2015）は、広く使われる計算ソフトウェアとして論文で公開 | 出典10 https://doi.org/10.1016/j.cpc.2021.108171 出典11 https://doi.org/10.1016/j.softx.2015.06.001 | 事実 |
| 分子シミュレーションの基礎の指針とチェックリスト（Braun ら）。不確かさとサンプリング品質の指針（Grossfield ら）。後者は、力場の不正確さのような系統的な誤差を扱わないと明記 | 出典32 https://doi.org/10.33011/livecoms.1.1.5957 出典33 https://doi.org/10.33011/livecoms.1.1.5067 | 事実 |
| 素材技術者の視点（対応表、力場の検証、温度で速めることの代償、初期構造の記録）。公開された力場が代表的な系で作られているという記述は推測 | 解説（本稿の見解） | |
| MDの力場・手順を定める国際規格は確認できなかった。時間刻み・力場の認証制度を定める規格も確認できなかった | 本稿の確認結果（確認できなかった） |
最終更新：2026年10月22日／出典は一次情報（原論文、ソフトウェアの公式文書、NIST・OpenKIM・CEN・ノーベル財団の公式ページ）に限定しています。数値例は本稿の計算で、特定の物質の実測値ではありません。Allen と Tildesley の教科書の個々の記述、C–H 伸縮の波数、簡易MDの文献値との比較、実際の素材での力場の精度、MDの力場や手順を定める国際規格については、確認できなかったため、断定していません。図2〜図8は説明用の概念図または本稿の計算を含む作図で、図1と図9はAI生成イメージ（描く場面の指示のみ記載）です。