応力鎖

応力鎖（おうりょくさ、英語: force chain）は、粒状体や粉体材料において、粒子間接触を介して伝達される力が空間的に不均一に集中し、鎖状・網状の構造を形成したものを指す。

離散要素法（Discrete Element Method、 DEM）では、接触力と粒子応力テンソルを計算することで、応力鎖ネットワークを定量的に抽出し、材料の強度・変形・流動性などの巨視的特性との関係を評価することができる。

概要

粒状体内部では、外力が一様に分布するのではなく、接触力が一部の粒子列に集中して伝達される。このような高応力の粒子列が「応力鎖」であり、それ以外の粒子群は比較的弱い接触力で全体を支持する「弱いネットワーク」として振る舞うと理解されている。

実験的には、フォトエラスティック円盤などを用いた可視化により、光弾性しまの強度が高い準一次元的な経路として観察され、数値的には DEM による接触力解析から同様の構造が確認されている。応力鎖の長さ分布や空間配置は、せん断帯の形成、局所的な座屈・破壊、流動開始（ジャミング／アンジャミング）などと密接に関連することが報告されている。

定義と理論

粒子応力とマクロ応力

連続体力学における応力テンソル（Cauchy 応力） $\sigma _{ij}$ は、単位面に作用する力ベクトルと面法線ベクトルの関係で定義される。DEM では、粒子同士の離散的な接触力からマクロな応力テンソルを構成するために、Love–Weber 型のマイクロ–マクロ関係式が広く用いられる。

粒子集団が体積 $V$ を占め、その中の全接触を $c=1,\dots ,N_{c}$ とすると、マクロ応力テンソルは

$\sigma _{ij}={\frac {1}{V}}\sum _{c=1}^{N_{c}}f_{i}^{(c)}\,\ell _{j}^{(c)}$

と表される。ここで $f_{i}^{(c)}$ は接触 $c$ に作用する接触力ベクトルの $i$ 成分、 $\ell _{j}^{(c)}$ は粒子中心から接触点への「ブランチベクトル」の $j$ 成分である。

個々の粒子 $p$ に対しても同様に、粒子近傍の体積を $V_{p}$ 、粒子 $p$ に属する接触集合を ${\mathcal {C}}_{p}$ とすると、粒子応力テンソル

$\sigma _{ij}^{(p)}={\frac {1}{V_{p}}}\sum _{c\in {\mathcal {C}}_{p}}f_{i}^{(c)}\,\ell _{j}^{(c)}$

を定義できる。このテンソルの最大固有値 $\lambda _{\max }^{(p)}$ はしばしば粒子が担う「粒子応力」や「荷重強度」を表すスカラー指標として用いられる。

応力鎖の定量的定義

応力鎖の定義にはさまざまな提案があるが、離散要素解析に基づく定量的な定義の一例は以下のように整理される。

応力鎖を構成する粒子は、粒子応力 $\lambda _{\max }^{(p)}$ が集団平均 $\langle \lambda _{\max }\rangle$ より大きい（あるいはしきい値 $\alpha \langle \lambda _{\max }\rangle$ を超える）粒子とする。
粒子中心を結ぶブランチベクトルと主応力方向（ $\lambda _{\max }^{(p)}$ に対応する固有ベクトル）とのなす角が小さい粒子対を「同一鎖上の隣接粒子」とみなす。
上記条件を満たす粒子が 3 個以上、準一次元的に連なった集合を応力鎖と定義する。

このような定義に基づけば、DEM 計算から得られた接触力と粒子配置に対し、(1) 粒子応力テンソルの計算、(2) 主応力方向の算出、(3) しきい値・角度条件による粒子グラフ上でのクラスタリング、といった手続きで応力鎖ネットワークを抽出できる。

数値シミュレーションにおける応力鎖

DEM モデルと接触力

DEM では、粒子間相互作用をバネ–ダッシュポットモデルや Hertz–Mindlin モデルなどの接触モデルで表現し、各タイムステップで粒子運動方程式を数値積分する。粒子間接触力は通常、法線方向成分 $\mathbf {F} _{\mathrm {n} }$ と接線方向成分 $\mathbf {F} _{\mathrm {t} }$ に分解され、摩擦や転がり抵抗などもモデル化される。

応力鎖解析では、これらの接触力ベクトルと粒子位置から、前述の Love–Weber 型関係式で粒子応力・マクロ応力を計算し、高応力粒子および強い接触のネットワークを抽出する。

LIGGGHTS / LAMMPS による粒子ダンプの出力

オープンソース DEM コードである LIGGGHTS や LAMMPS では、dump コマンドにより粒子データをファイル（dump ファイル）に出力できる。粒子位置や速度、力、半径などを含むダンプの例を以下に示す。

# 粒子データのダンプ例（LIGGGHTS / LAMMPS）
dump    dp all custom 400 dump.p &
        id type x y z ix iy iz vx vy vz fx fy fz omegax omegay omegaz radius

このような粒子ダンプは、VTK 形式に変換して ParaView などで可視化することで、粒子ごとの応力指標を色付けした散布図や、接触力の大きさによるフィルタリングなどに利用できる。

接触局所量（pair/local, dump local）の出力

応力鎖ネットワークを明示的に描画するには、粒子–粒子または粒子–壁の「接触ごとの局所量」を出力する必要がある。LIGGGHTS では compute pair/gran/local（および compute wall/gran/local）コマンドにより、個々の接触について位置・ID・接触力などを配列として計算し、dump local でファイルに書き出すことができる。

粒子–粒子接触に対し、接触している 2 粒子の位置、それぞれの ID、周期境界フラグ、全接触力、法線・接線成分を出力する設定の例を示す。

# 接触局所量の計算（粒子–粒子）
compute cp all pair/gran/local id pos force force_normal force_tangential

# 接触ごとのローカルダンプ（18 列）
dump    fc all local 400 post/forcechain-*.dump  &
        c_cp[1]  c_cp[2]  c_cp[3]  c_cp[4]  c_cp[5]  c_cp[6]   &  # x1 y1 z1 x2 y2 z2
        c_cp[7]  c_cp[8]  c_cp[9]                                 # id1 id2 periodic
        c_cp[10] c_cp[11] c_cp[12]                                # fx  fy  fz
        c_cp[13] c_cp[14] c_cp[15]                                # Fnx Fny Fnz
        c_cp[16] c_cp[17] c_cp[18]                                # Ftx Fty Ftz

このような dump ファイルは、行ごとに 1 つの接触を表し、列として

粒子 1 の座標 $(x_{1},y_{1},z_{1})$
粒子 2 の座標 $(x_{2},y_{2},z_{2})$
粒子 1, 粒子 2 の ID と周期境界フラグ
全接触力 $\mathbf {F} =(F_{x},F_{y},F_{z})$
法線成分 $\mathbf {F} _{\mathrm {n} }=(F_{\mathrm {n} x},F_{\mathrm {n} y},F_{\mathrm {n} z})$
接線成分 $\mathbf {F} _{\mathrm {t} }=(F_{\mathrm {t} x},F_{\mathrm {t} y},F_{\mathrm {t} z})$

などを含む。粒子–壁接触についても compute wall/gran/local と dump local を組み合わせることで、同様の形式で出力できる。

dump ファイルの前処理

ENTRIES ヘッダーと列名の付与

LIGGGHTS の dump local は、接触ローカル量を

ITEM: ENTRIES id fx fy fz ...

のように、LIGGGHTS 固有の列順で出力する。しかし、後処理ツールによっては、「1 列目は x1、2 列目は y1 …」といった物理量の名前付きヘッダーを期待する場合がある。

この問題に対して、Python スクリプトを用いてヘッダーを書き換える方法がよく用いられる。例えば、GitHub 上で公開されている liggghts-entries-header-replacer は、LIGGGHTS dump ファイル中の

ITEM: ENTRIES ...

行を任意の 1 行に一括置換するコマンドラインツールであり、標準ライブラリのみで動作する。使い方の概略は次のとおりである。

# リポジトリを取得
git clone https://github.com/opendemjapan/liggghts-entries-header-replacer
cd liggghts-entries-header-replacer

# ヘッダーの一括置換（別名保存）
python replace_entries_header.py \
    -H "ITEM: ENTRIES x1 y1 z1 x2 y2 z2 id1 id2 periodic fx fy fz Fnx Fny Fnz Ftx Fty Ftz" \
    post/forcechain-*.dump

この例では、post/forcechain-*.dump にマッチするファイルが、拡張子の手前に -renamed を付加した新しいファイル名で出力され、ITEM: ENTRIES ... 行だけが指定したヘッダーに置き換えられる。--inplace オプションを用いれば、安全な一時ファイル置換を伴う上書き処理も可能である。

dump から VTK への変換

Python スクリプトによる変換（liggghts-lammps-dump2vtk）

DEM のダンプを可視化するには、VTK 形式への変換が一般的である。従来、LAMMPS / LIGGGHTS 用の後処理プログラムである Pizza.py / LPP は Python 2 系で実装されており、Python 3 環境では利用しづらいという問題があった。

この制約を解消するため、粉体粒子（ふんたいつぶこ、Open DEM Japan）により、LIGGGHTS / LAMMPS の粒子ダンプ（dump.* や dump.*.gz）を VTK Legacy 形式の .vtk ファイルに一括変換する Python 3 スクリプト liggghts-lammps-dump2vtk が作成・公開されている。このスクリプトは、旧 Pizza.py（lpp.py / vtk.py / dump.py）に相当する機能を Python 3 向けに統合し、

単一ファイルのスクリプトとして配布されていること、
VTK Legacy の ASCII / BINARY 出力や、多数タイムステップのチャンク分割・並列処理に対応すること、

などを特徴とする。出力される .vtk ファイルには、粒子応力や接触数といったスカラー量、粒子速度や力といったベクトル量を点データとして格納でき、ParaView や Mayavi などの可視化ソフトウェアで直接読み込める。

Windows GUI / CLI ツール dump2vtk.exe（D2V）

上記の Python 3 スクリプトをベースにした Windows 用スタンドアロンツールとして、dump2vtk.exe（通称 D2V）が配布されている^[1]。開発者は前節と同じく粉体粒子であり、Python 3 版 LPP/Pizza を誰でも簡単に利用できるよう、PySide6 ベースの GUI と同等機能を備えた CLI を備える形で提供している。D2V は LIGGGHTS / LAMMPS の dump を VTK（レガシー形式、ASCII / BINARY）へ高速一括変換するツールであり、Windows 用バイナリ（インストーラ形式）として配布されている。

GUI は 3 つのタブから構成される。

Particles – 粒子ダンプを VTK に変換するタブ。1 つの dump ファイルをドラッグ＆ドロップすると、同じプレフィクスを共有する連番ファイルを自動収集し、並列処理で一括変換する。
Header Replace – ITEM: ENTRIES ... 行をバッチ置換するタブ。前節のヘッダー置換スクリプトと同等の機能を GUI から提供する。
Force chain – pair/local 形式の dump からフォースネットワーク（force chain）を VTK として出力するタブ。Louvain 法によるコミュニティ検出により、力のネットワークをクラスタごとに色分けするテスト的な解析機能も備える。

特に Force chain タブは、GUI 上では「応力」や「応力鎖」の可視化として紹介されることがあるが、実際には pair/local dump に記録された接触ごとの力を基に、強い接触を抽出して線分メッシュとして描画する「力のネットワーク（force chain）」の可視化が主目的である。前節で述べた Love–Weber 型の粒子応力テンソルを通じて Cauchy 応力と対応付けることはできるものの、D2V 自体は連続体力学における応力場全体を厳密に再構成する解析コードではなく、「どの経路に大きな力が流れているか」を可視化するツールとして利用されることが多い。

CLI（CUI）では、GUI の各タブに対応する機能がサブコマンドとして提供されており、バッチ処理や自動実行に利用できる。代表的な用途としては、

LPP/Pizza 互換の粒子 dump → VTK 変換、
ITEM: ENTRIES ... ヘッダーの一括置換、
pair/local dump からの force chain VTK 出力（Louvain 法によるコミュニティ集約量の出力を含む）、

などが挙げられる。force 系サブコマンドは、前処理済みの pair/local dump（例えば *-renamed.dump）を入力として、タイムステップごとに接触ネットワークの VTK を出力する。これにより、粒子中心を結ぶ線分として応力鎖（力のネットワーク）を可視化したり、接触力の大きさに応じて線の太さや色を変えたりすることが容易になる。

可視化と解析

VTK 形式に変換された粒子・接触データは、ParaView などの可視化ツールで以下のような解析に利用される。

粒子を点として描画し、粒子応力や接触数、主応力方向などのスカラー／ベクトル場を色付け・矢印表示する。
接触ダンプから生成された線分メッシュを利用し、接触力の大きさに応じた太さ・色で「力のネットワーク」を描画する。
一定しきい値以上の接触力だけを表示して「強いネットワーク（strong force network）」を抽出し、応力鎖の空間分布や長さ分布を評価する。
Louvain 法などのグラフクラスタリングにより、応力鎖ネットワークをコミュニティ別に色分けし、局所構造とマクロ応答の関係を調べる。

このような可視化と統計解析により、接触力・応力鎖・マクロな荷重–変位応答の三者の相関を体系的に調査することが可能となる。

脚注

^ D2V

外部リンク

LIGGGHTS 公式サイト（英語） – https://www.cfdem.com/liggghts-open-source-discrete-element-method-particle-simulation-code
Open DEM Japan – Python による LIGGGHTS ワークフローと dump2vtk 解説
GitHub – opendemjapan/liggghts-lammps-dump2vtk
GitHub – opendemjapan/liggghts-entries-header-replacer
Open DEM Japan – 「Python による LIGGGHTS® ワークフロー」（dump2vtk GUI 配布ページ） – https://opendemjapan.parallel.jp/python.html

[1] D2V

[1]

概要