個別要素法による粒状体群のせん断シミュレーションにおける摩擦処理

70
法政大学情報メディア教育研究センター研究報告 Vol.29
2015 年
個別要素法による粒状体群のせん断シミュレーションにおける摩擦処理
On the Frictional Contact in Shear Simulation of Particles by DEM
板谷知洋１）大西泰史２）吉田長行２）
Tomohiro Itadani, Yasushi Onishi, Nagayuki Yoshida
１）
法政大学大学院デザイン工学研究科建築学専攻
法政大学デザイン工学部建築学科
２）
In the Distinct Element Method (DEM), Friction is introduced in the tangential direction between the two
particles. The way of its modeling could greatly influence on the analytical result. Therefore the ordinal calculation
scheme have to be reconsidered in the tangential direction of the contact surface. In this paper, the numerical
procedure taking account of the stick-slip phenomenon is introduced to deal with the friction among particles more
precisely. The analytical result shows that smoothness of force and displacement relation can be improved.
Furthermore it is confirmed that the random size on particles is relatively more stable than the uniform size.
Keywords : DEM, the direct shear test, tangential force, Stick-Slip phenomenon
1. はじめに
2. 解析手法
不連続体を扱う手法の一つに Cundall.P.A[1],[2]が
提唱した個別要素法（Distinct Element Method：以降
DEM と略称する）が挙げられる。この手法の特徴は
粒状体の個々の粒子に働くミクロな力の相互作用を
考慮することで材料のマクロな力学挙動を再現でき
るところにあり、連続体力学とは別の枠組みとなる
手法である。この手法を用いれば、本来大がかりな
機材や時間を有する土質試験をコンピュータ上でシ
ミュレーションすることが可能である。また、その
際に粒状体に生じる進行性破壊、ダイレイタンシー
現象、せん断帯形成過程の視覚的な追跡が容易に行
える。
粒状体群に働く接線方向力には摩擦が関与し、解
析結果に大きな影響を与える。そこで本研究では、
ベルトコンベアー上にバネを介して置いた物体に繰
り返し見られる Stick-Slip 現象を粒子間の相対回転
に取り入れ、計算手順を定式化し、せん断シミュレ
ーションを通じて新たな手法の有用性を確認する。
原稿受付 2015 年 3 月 9 日
発行
2015 年 4 月 1 日
Copyright © 2015 Hosei University
2.1 DEM
DEM の計算手法は伯野の文献[3]および文献[4]に
詳しい。本論では粒子間の接線方向における摩擦の
取り扱いをより厳密化するため Stick-Slip 現象を考
慮した計算手法を提案し、これを Stick-Slip 型と呼称
する。
2.2 Stick-Slip 型計算手法
Stick-Slip 型では法線方向力の計算手順は変更せ
ず、接線方向力及びこの力によって起こるトルクの
評価法を再検討する。
先ず、2 つの粒子 i、j が接触状態のとき(1)Stick 状
態か(2)Slip 状態かの判定を行い、それぞれの計算に
入る。以下に計算法を示す。
(1) Stick 状態
ある時刻の間、粒子 i が粒子 j に対して滑らない
で回転する場合、粒子間の接線方向相対変位の回転
71
成分は次の条件を満たす必要がある。
u t  ri  i t  rj   j t  0
(1)
しかし、この間滑らずに粒子の面上を転がるために
は接線方向力を作用させる必要がある。そこで、あ
る時刻に生じる相対変位を並進による相対速度より
求め、接線方向力を次式で与える。
Tij   ks uT t  t  cs uT t
t
(2)
この力に対し生じる回転量についてはモーメントの
式を解かず、得られた増分並進変位から求め次式と
する。
  
i t

 uT t
(3)
ri
(2) Slip 状態
粒子 i が粒子 j に対して滑って回転する場合、接
線方向並進変位はないと考える。しかし滑る以上は
並進せずに回転しているのでトルクの評価が必要で
ある。トルクを発生させる接線方向力を滑る瞬間の
値として、次式で与える。




 M i t  sign us t     Nij t  cs  us t  ri
(4)
3. 一面せん断試験シミュレーション
3.1 モデル作成
初期配置を作成するために一時的に設置した箱に
自然落下により充填する。その様子を Fig.1 に示す。
図.2 締固め
Fig.2 Compaction
3.2 解析内容
本解析に用いたパラメータを Table.1、解析条件を
Table.2 と Table.3 に示す。
表 1 物性値
Table.1 Analytical data
法線バネ定数
kn  1.0×106 [N/cm]
接線バネ定数
ks  2.5×105 [N/cm]
内部摩擦角
  27 [deg]
  0.51[]
静止摩擦係数
粒子密度
減衰定数(落下時)
  2.65×103 [kg / cm3 ]
hn , hs  1.0[]
減衰定数
hn , hs  0.215[]
時間増分
せん断速度
t  1.0×106 [sec]
Vx  1.0 [mm/sec]
Stick-Slip 型を用いる場合、Slip 状態時のトルクの評
価を式(4)で与える。その際に滑る瞬間の値を動摩擦
係数で計算する。本解析では動摩擦係数を静止摩擦
係数の 40%と 60%の 2 つの値を設定する。
1  0.4
2  0.6
Pattern1
Pattern2
Pattern3
Pattern4
図.1 充填
Fig.1 Packing
(5)
表 2 解析条件 1
Table.2 Analytical condition1
粒径[cm] 粒子数[個] 動摩擦係数
1
0.2
1000
1
0.4
1000
1
0.2~0.4
1000
1
0.2~0.8
1000
脊黒ら[5]が導入した周期境界、固定粒子境界を用
いて、粒子が安定状態に移行するまで時間を置く。
その様子を Fig.2 に示す。
Copyright © 2015 Hosei University
法政大学情報メディア教育研究センター研究報告 Vol.29
72
Table.3 Analytical condition2
粒径[cm] 粒子数[個] 動摩擦係数
 2
0.2
4000
 2
0.4
4000
 2
0.2~0.4
4000
 2
0.2~0.8
4000
Pattern5
Pattern6
Pattern7
Pattern8
Shear stress[kN/㎡]
表 3 解析条件 2
400
300
200
100
0
0
せん断応力の評価は底板に作用する水平方向の力
の総和を解析領域幅で除して算出する。
SFx   f x
(6)
4. 解析結果
4.1 せん断応力、ダイレイタンシー
Table.2 と Table.3 の条件で解析を行い、せん断応
力、
ダイレイタンシー曲線を比較検討する。
各 Pattern
拘束圧は 100[kPa]、200[kPa]、300[kPa]、400[kPa]の
4 段階で検証し、グラフには順に赤、緑、紫、青で
色分けしている。以下に Table.2 と Table.3 の結果を
示す。
2
0.2
Dilatancy[cm]
 f x ：せん断板に作用する水平方向力
0.5
1
1.5
Displacement[cm]
0.15
0.1
0.05
0
-0.05
0
0.5
1
1.5
2
Displacement[cm]
図.4 粒径 0.4cm・粒子数 1000 個
Fig.4 Particle size 0.4[cm] and Particle number 1000
Shear stress[kN/㎡]
250
Shear stress[kN/㎡]
500
400
300
200
100
200
150
100
50
0
-50
0
-100 0
0.5
1
1.5
Displacement[cm]
0
0.5
1
1.5
Displacement[cm]
2
2
0.15
Dilatancy[cm]
Dilatancy[cm]
0.2
0.15
0.1
0.05
0.1
0.05
0
0
0
-0.05
-0.1
0
0.5
1
1.5
2
Displacement[cm]
図.3 粒径 0.2cm・粒子数 1000 個
Fig.3 Particle size 0.2[cm] and Particle number 1000
Copyright © 2015 Hosei University
-0.05
0.5
1
1.5
Displacement[cm]
図.5 粒径 0.2~0.4cm・粒子数 1000 個
Fig.5 Particle size 0.2~0.4[cm] and
Particle number 1000
法政大学情報メディア教育研究センター研究報告 Vol.29
2
73
Shear stress[kN/㎡]
200
150
100
50
0
-50
Dilatancy[cm]
1500
0
0.6
0.5
0.4
0.3
0.2
0.1
0
-0.1 0
0.5
1
1.5
Displacement[cm]
2
0.5
1
1.5
Displacement[cm]
0.5
1.0
1.5
2.0
Displacementr[cm]
0.8
0.6
0.4
0.2
0
300
200
100
0
0
0.5
1
1.5
0.5
1
1.5
Displacement[cm]
2
図.8 粒径 0.4cm・粒子数 4000 個
Fig.8 Particle size 0.4[cm] and Particle number 4000
Shear stress[kN/㎡]
Shear stress[kN/㎡]
0.0
2
250
200
150
100
50
0
2
0
0.5
1
1.5
Displacement[cm]
2
0
0.5
1
1.5
Displacement[cm]
2
Displacement[cm]
0.1
Dilatancy[cm]
0.2
Dilatancy[cm]
0
0
400
0.15
0.1
0.05
0
-0.05
500
1
図.6 粒径 0.2~0.8cm・粒子数 1000 個
Fig.6 Particle size 0.2~0.8[cm] and
Particle number 1000
-100
1000
-500
Dilatancy[cm]
Shear stress[kN/㎡]
250
0.08
0.06
0.04
0.02
0
0
0.5
1
1.5
2
Displacement[cm]
図.7 粒径 0.2cm・粒子数 4000 個
Fig.7 Particle size 0.2[cm] and Particle number 4000
Copyright © 2015 Hosei University
図.9 粒径 0.2~0.4cm・粒子数 4000 個
Fig.9 Particle size 0.2~0.4[cm] and
Particle number 4000
法政大学情報メディア教育研究センター研究報告 Vol.29
74
Table.3
600
200
100
0
0
0.5
-100
1
1.5
2
Shear stress[kN/㎡]
Shear stress[kN/㎡]
300
500
400
300
200
100
Displacement[cm]
0
0
100
200
300
400
Shear stress[kN/㎡]
Normal stress[kN/㎡]
5
4
図.12 せん断抵抗角
Fig.12 Angle of shear resistance
3
2
5. 考察・結論
1
0
0
0.5
1
1.5
Displacement[cm]
2
図.10 粒径 0.2~0.8cm・粒子数 4000 個
Fig.10 Particle size 0.2~0.8[cm] and
Particle number 4000
4.2 せん断強度、せん断抵抗角
せん断応力図から各拘束圧ごとにピークをプロッ
トし、垂直応力－せん断応力関係の図からせん断抵
抗角を算出する。Pattern 順に赤、緑、紫、青に色分
けしている。以下に Table.2 と Table.3 の結果を示す。
Table.2
Shear stress[kN/㎡]
400
粒径 0.2～0.4cm では比較的安定した結果が見られ
る。粒径範囲を広くした 0.2～0.8cm では粒子の大小
の差が大きいため粒子同士のかみ合いが頻繁に起こ
り，せん断応力の急激な減少を引き起こしている。
また、粒子数 1000 個と 4000 個の結果を比較すると
目立った差はないため，解析時間短縮のためより少
ない粒子数での解析においても十分な結果が得られ
ることが見込める。
6. 結論
Stick-Slip 型を導入した場合、解析時間の増加は見
られず、安定した結果も得ることができた。また、
精度を保証するための粒子数の下限値を推定した。
参考文献
300
200
100
0
0
100
200
300
Normal stress[kN/㎡]
図.11 せん断抵抗角
Fig.11 Angle of shear resistance
Copyright © 2015 Hosei University
400
[1] Cundall,P,A.：A computer model for simulating
progressive, large-scale movements in blocky
rock system. ISRM，Nancy，France. Proc. ，2，
129-136，1971
[2] Cundall,P.A. STRACK, O.D.L. ： A discrete
numerical model for granular assemblies,
Géotechnique，29，1，47-65，1979
[3] 伯野元彦：
「破壊のシミュレーション」―拡張
個別要素法で破壊を負う―，森北出版，1997
[4] 粉体工学会編：粉体シミュレーション入門―
コンピュータで粉体技術を創造する，産業図
書，1998
[5] 脊黒隆大、板谷知洋、吉田長行：土質試験シ
法政大学情報メディア教育研究センター研究報告 Vol.29
75
ミュレーションにおける個別要素法の特性検
討、法政大学情報メディア教育センター研究
報告、Vol.28、2014
Copyright © 2015 Hosei University
法政大学情報メディア教育研究センター研究報告 Vol.29

Download Report