剛塑性体における平面ひずみ問題の例に，よく知られているブロックの圧縮とポンチの押し込み問題がある．一般に，金属板より圧縮板の幅が大きいときその問題はブロックの圧縮問題として解析され，その逆の場合はポンチの押し込み問題として解析される．
$\hspace{0em}$図7.1 ポンチの押し込み問題
本節では
図7.1 に示すポンチの押し込み問題に対する数値計算例を示す．ポンチの押し込み問題は解析的に考えると，ブロックの圧縮問題と異なり，一様な圧縮応力場と剛塑性域の境界上で速度の不連続がゆるされず，取り扱いが面倒になる．しかし，数値計算では境界条件以外の差はほとんど無く，構造が複雑になればなるほど解析的手法に比べ威力を発揮する．
本書で用いた寸法及び要素分割を
図7.2 に示す．この例題では $2h=5{\rm cm}, 2b=2.5{\rm cm}, h/b=2$ としている．ただし，
図7.1 に示すように左右対称であるため，全体の左側半分を解析モデルとして設定している．節点数が15節点，領域数（要素数）が18要素で，総自由度数は108自由度である．
$\hspace{5em}$図7.2 要素分割例
拘束条件としては，領域番号⑫の辺7-14の上下方向と，領域⑫の14-11，領域⑭の11-12，領域⑯の12-13，領域⑱の13-15の左右方向を拘束する．また．領域⑱における辺10-15に下向きの強制変位として
\[
\hat{\delta} = -1.0 \times 10^{-3} {\rm m}
\]
を与えて変位を拘束する．
破壊条件としてはMohr-Coulombの条件を用い，材料定数を
表7.1 に示すよう仮定した．ただし，内部摩擦角を0としているため，実際にはTrescaの条件を適用したことになる．また，応力状態しては平面ひずみ状態を仮定している．
表7.1 材料定数
材料種類
値
弾性係数 (E)
200 ${\rm GPa}$ ポアソン比 ($\nu$)
0.3 せん断強さ (c)
3 ${\rm MPa}$ 内部摩擦角 ($\phi$)
$0^{\circ}$
なお，本例題では自重を無視している．上記の条件を基に作成された本書のサンプルプログラム用のデータは
図7.3 に示されている．
----+----1----+----2-----
sample data
5
15
1 0. 0.
2 0. 0.025
3 0. 0.05
4 0.0125 0.0125
5 0.0125 0.025
6 0.0125 0.0375
7 0.025 0.
8 0.0325 0.016
9 0.0325 0.034
10 0.025 0.05
11 0.0375 0.0125
12 0.0375 0.025
13 0.0375 0.0375
14 0.0375 0.
15 0.0375 0.05
18 0
1 1 2 1 4
2 1 2 4 5
3 1 2 5 6
4 1 6 3 2
5 1 1 7 4
6 1 4 7 8
7 1 5 4 8
8 1 8 9 5
9 1 6 5 9
10 1 9 10 6
11 1 10 3 6
12 1 7 14 11
13 1 8 7 11
14 1 8 11 12
15 1 9 8 12
16 1 9 12 13
17 1 10 9 13
18 1 10 13 15
1
1 200000000. 0.3 0. 3000. 0. 1.
0
6
1 12 14 7 2 0.
2 12 11 14 1 0.
3 14 12 11 1 0.
4 16 13 12 1 0.
5 18 15 13 1 0.
6 18 10 15 2 -1.
----+----1----+----2-----
図7.3 入力データ例
図7.4 ～
図7.11 は本書のサンプルプログラムによる入力情報の出力結果である．
図7.4 は基本データに関する情報を出力したもので，要素タイプは平面ひずみであり，プログラム内部で作成されたペナルティ辺の数は22である．また，メカニズム・チェック用の最大変位は1000としており，最大バンド幅は48である．
TITLE : sample data
INFORMATION :
ELEMET TYPE ------------------> PLANE STRAIN
NUMBER OF NODE ---------------> 15
NUMBER OF ELEMENTS -----------> 18
NUMBER OF MATERIALS ----------> 1
NUMBER OF SPRINGS ------------> 22
NUMBER OF LOADS --------------> 0
NUMBER OF SUPPORT D.O.F. -----> 6
MAXIMUM BAND WIDTH -----------> 48
NUMBER OF ITARATIONS ---------> 5
MAXIMUM DISPLACEMENT ---------> 0.1E+04
図7.4 入力情報の出力結果
図7.5 は入力した座標値を出力したものである．
COORDINATE DATA :
NODE NO. X-COORDINATE Y-COORDINATE
1 0.00000E+00 0.00000E+00
2 0.00000E+00 0.25000E-01
3 0.00000E+00 0.50000E-01
4 0.12500E-01 0.12500E-01
5 0.12500E-01 0.25000E-01
6 0.12500E-01 0.37500E-01
7 0.25000E-01 0.00000E+00
8 0.32500E-01 0.16000E-01
9 0.32500E-01 0.34000E-01
10 0.25000E-01 0.50000E-01
11 0.37500E-01 0.12500E-01
12 0.37500E-01 0.25000E-01
13 0.37500E-01 0.37500E-01
14 0.37500E-01 0.00000E+00
15 0.37500E-01 0.50000E-01
図7.5 座標データの出力結果
図7.6 は小領域に関係するデータの出力結果である．出力内容は左から，
要素番号
要素構成節点番号(1)(2)(3)
材料種類番号(MAT.NO.)
要素面積(AREA)
x方向図心座標(XG)
y方向図心座標(YG)
の順に表示されている．要素構成節点番号(1)(2)(3)は先に示したように，三角形の頂点に関する節点番号を左回りに読み込んだデータである．
ELEMENT DATA :
ELEM.NO. (1) (2) (3) MAT.NO. AREA XG YG
1 2 1 4 1 0.15625E-03 0.41667E-02 0.12500E-01
2 2 4 5 1 0.78125E-04 0.83333E-02 0.20833E-01
3 2 5 6 1 0.78125E-04 0.83333E-02 0.29167E-01
4 6 3 2 1 0.15625E-03 0.41667E-02 0.37500E-01
5 1 7 4 1 0.15625E-03 0.12500E-01 0.41667E-02
6 4 7 8 1 0.14688E-03 0.23333E-01 0.95000E-02
7 5 4 8 1 0.12500E-03 0.19167E-01 0.17833E-01
8 8 9 5 1 0.18000E-03 0.25833E-01 0.25000E-01
9 6 5 9 1 0.12500E-03 0.19167E-01 0.32167E-01
10 9 10 6 1 0.14688E-03 0.23333E-01 0.40500E-01
11 10 3 6 1 0.15625E-03 0.12500E-01 0.45833E-01
12 7 14 11 1 0.78125E-04 0.33333E-01 0.41667E-02
13 8 7 11 1 0.53125E-04 0.31667E-01 0.95000E-02
14 8 11 12 1 0.31250E-04 0.35833E-01 0.17833E-01
15 9 8 12 1 0.45000E-04 0.34167E-01 0.25000E-01
16 9 12 13 1 0.31250E-04 0.35833E-01 0.32167E-01
17 10 9 13 1 0.53125E-04 0.31667E-01 0.40500E-01
18 10 13 15 1 0.78125E-04 0.33333E-01 0.45833E-01
図7.6 要素データの出力結果
プログラム内部で作成したペナルティ関係のデータ出力は，
図7.7 に示す記号にしたがい，左から以下の順に表示されている．このときの出力結果は，
図7.8 に示されている．
$\hspace{0em}$図7.7 ペナルティデータに関する記号
ばね番号
ばね構成要素番号(E-1)(E-2)
ばね構成節点番号(N-1)(N-2)
SPRING DATA :
SPRG.NO. E-1 E-2 N-1 N-2 L H1 H2
1 1 2 4 2 0.17678E-01 0.58926E-02 0.29463E-02
2 1 5 1 4 0.17678E-01 0.58926E-02 0.58926E-02
3 2 3 5 2 0.12500E-01 0.41667E-02 0.41667E-02
4 2 7 4 5 0.12500E-01 0.41667E-02 0.66667E-02
5 3 4 6 2 0.17678E-01 0.29463E-02 0.58926E-02
6 3 9 5 6 0.12500E-01 0.41667E-02 0.66667E-02
7 4 11 6 3 0.17678E-01 0.58926E-02 0.58926E-02
8 5 6 7 4 0.17678E-01 0.58926E-02 0.55390E-02
9 6 7 8 4 0.20304E-01 0.48225E-02 0.41043E-02
10 6 13 7 8 0.17671E-01 0.55412E-02 0.20043E-02
11 7 8 8 5 0.21932E-01 0.37997E-02 0.54715E-02
12 8 9 9 5 0.21932E-01 0.54715E-02 0.37997E-02
13 8 15 8 9 0.18000E-01 0.66667E-02 0.16667E-02
14 9 10 9 6 0.20304E-01 0.41043E-02 0.48225E-02
15 10 11 10 6 0.17678E-01 0.55390E-02 0.58926E-02
16 10 17 9 10 0.17671E-01 0.55412E-02 0.20043E-02
17 12 13 11 7 0.17678E-01 0.29463E-02 0.20035E-02
18 13 14 11 8 0.61033E-02 0.58029E-02 0.34135E-02
19 14 15 12 8 0.10296E-01 0.20235E-02 0.29139E-02
20 15 16 12 9 0.10296E-01 0.29139E-02 0.20235E-02
21 16 17 13 9 0.61033E-02 0.34135E-02 0.58029E-02
22 17 18 13 10 0.17678E-01 0.20035E-02 0.29463E-02
図7.8 ペナルティデータの出力結果
図7.9 は材料データであり，対応している出力内容は左から以下の順で表示されている．本例題では平面ひずみを仮定しているため，解析板厚を1としている．
材料種類番号
弾性係数(E)
ポアソン比(P)
物体力(BODY FORCE)
せん断強さ(C)
内部摩擦角(PHI)
解析板厚(t)
MATERIAL DATA :
MAT.NO. E P BODY FORCE C PHI t
1 0.20000E+09 0.30000E+00 0.00000E+00 0.30000E+04 0.00000E+00 0.10000E+01
図7.9 材料データの出力結果
図7.10 は載荷データで載荷要素番号の他，載荷値が以下のように出力される．本サンプルプログラムでは，領域の境界辺に分布荷重が作用する場合のみを取り扱っているため，領域番号に続いて，荷重が作用する境界辺を構成する節点番号が続いて表示される．
図7.10 は荷重の出力例である．ただし，図の結果は荷重を作用させた場合のみ出力されるもので，本章で取り上げている例題は，強制変位のみ作用する，いわゆる変制御による解析例であるため，
図7.10 の出力は表示されない．
荷重番号
載荷領域番号
節点1
節点2
荷重値(X)
荷重値(Y)
LOAD DATA :
NO. ELEM.NO. N-1 N-2 X Y
1 4 6 3 0.100E+01 0.000E+00
図7.10 荷重データ出力結果
図7.11 は拘束条件で，対応している出力内容は左から以下の順で表示されている．本サンプル・プログラムでは，領域の境界辺の拘束のみ取り扱っているため，領域番号に続いて，拘束する領域の辺を構成する節点番号がつづいて表示される．拘束方向としてはは $x$ 方向が1，$y$ 方向が2と表示される．また，移動しない支点のような状態の場合は，強制変位が0と表示される．
拘束番号
拘束領域番号
節点1
節点2
拘束方向
強制変位
SUPPORT D.O.F.:
NO. ELEM.NO. N-1 N-2 DIREC. disp.
1 12 14 7 2 0.00000E+00
2 12 11 14 1 0.00000E+00
3 14 12 11 1 0.00000E+00
4 16 13 12 1 0.00000E+00
5 18 15 13 1 0.00000E+00
6 18 10 15 2 -0.10000E+01
図7.11 支点データ出力結果
本節では
7.1節 で説明したデータを基に，本サンプルプログラムにより解析した結果を示す．
HPMを用いた離散化極限解析における特徴の1つはすべり線が得られることである．従来のFEMでは，降伏した要素を塗りつぶしたりすることによって塑性化領域を表現していたが，すべり線の場合は小領域（要素）自身ではなく，小領域（要素）境界辺の塑性化であるため，FEMとは表現方法が異なる．いろいろな表現方法が考えられるが，本書では太線によりすべり線を表現している．
$\hspace{2em}$図7.12 すべり線の発達状況
$\hspace{4em}$図7.13 変位モードとすべり線
図7.12 は本例題に関するすべり線の発達状況を表したものである．図中の数字はすべりの発生した荷重増分ステップ番号を示している．本例題では，対称条件を用いているため，実際の解析結果は右側の網掛け部分を除く左側半分のみである．破壊は，載荷板近傍から発生し，中央部へと進行している．
HPMでは，領域内の剛性を評価しており，領域の変形を考慮することができる．
図7.13 は，変形状況とあわせてすべり線を描画した図である．ただし，$r_{\rm min}=1$ の場合と他の結果では変位スケールが異なっている．
また，HPMでは，領域内において，テンソル量としてのひずみや応力も求められる．この結果を利用すれば，
図7.14 のようにvonMisesの応力図を描いたりすることもできる．図はメカニズムを形成した，$r_{\rm min}=0.10535$ におけるvonMisesの応力図である．領域分割が粗いため，滑らかな分布ではないが，細分化することで，精度上げることができる．
$\hspace{1em}$図7.14 vonMiseの応力（ $r_{\rm min}=0.10535$ ）
計算された強制変位 と載荷荷重の関係は
表7.2 に整理されている．$p$ は強制変位を与えた辺での反力を載荷幅1.25cmで割つた単位面積当りの値である．また，表の最後に示した $p/2c$ は単位面積当りの載荷荷重をせん断強さにより無次元化した値である．
このときの極限荷重は解析的に
\[
{\rm (7.1)}
\frac{P_m}{2c} = 1.22
\]
表7.2 荷重増分率と載荷荷重
荷重増分ステップ
$\delta(\times 10^{-3} {\rm m})$
p (${\rm MPa}$)
p/2c
1
0.07935
6.415
1.07 2
0.08507
6.672
1.11 3
0.09192
6.966
1.16 4
0.10535
7.410
1.23
となる．ここで $p_m$ が載荷面上における単位面積当りの極限荷重である．粗い領域分割にもかかわらず，1%以下の誤差に収まっている．
一方，
表7.2 を基に作成した荷重－変位曲線を
図7.15 に示す．縦軸が $p/2c$ の値，横軸が強制変位量の値を示している．
$\hspace{2em}$図7.15 載荷板における荷重と変位の関係
最後に出力結果に対する説明を行う．荷重増分計算の場合，繰り返し計算が必要であり，本例題では4回の操り返し計算により崩壊荷重が得られている．その後，与えた強制変位量まで変形するため，出力は5ステップ分の量があるが，ここでは，メカニズム形成時のステップ4のみ取り上げて説明を行う．
図7.16 は各計算ステップにおける荷重増分率の値であり，出力は左から以下のような内容となっている．例題の出力では4ステップ目の荷重増分率が0.013431で，今回迄の累積荷重増分率が0.10535となっている．
STEP NO.= 荷重増分ステップ番号
STEP = 今回の荷重増分率
RMIN = 累積荷重増分率
OUTPUT DATA :
STEP NO.= 4 STEP = 0.13431E-01 RMIN = 0.10535E+00
図7.16 荷重増分率の出力結果
図7.17 は変位に関する出力結果である．結果は，
の順に表示される．単位は入力された単位系を採用している．
< DISPLACEMENT >
ELEM.NO. U V T
1 -0.21282E-06 -0.50357E-06 0.11179E-04
2 -0.27088E-06 -0.48434E-06 0.11977E-05
3 -0.27088E-06 -0.56921E-06 -0.11977E-05
4 -0.21282E-06 -0.54997E-06 -0.11179E-04
5 -0.13478E-06 -0.36212E-06 0.15160E-04
6 -0.12993E-06 -0.28250E-06 0.40214E-05
7 -0.18381E-06 -0.42394E-06 0.34359E-05
8 -0.13088E-06 -0.52677E-06 0.14014E-13
9 -0.18381E-06 -0.62960E-06 -0.34359E-05
10 -0.12993E-06 -0.77105E-06 -0.40214E-05
11 -0.13478E-06 -0.69142E-06 -0.15160E-04
12 -0.83619E-08 -0.77083E-07 0.34148E-08
13 -0.50783E-07 -0.21204E-06 0.40570E-05
14 -0.22662E-07 -0.37800E-06 0.97042E-06
15 -0.45317E-07 -0.52677E-06 -0.26338E-13
16 -0.22662E-07 -0.67554E-06 -0.97042E-06
17 -0.50783E-07 -0.84150E-06 -0.40570E-05
18 -0.83619E-08 -0.97646E-06 -0.34148E-08
図7.17 変位の出力結果
次に，領域内の応力とひずみが領域数だけ
図7.18 のように出力される．ただし，出力内容は，下記のように応力，ひずみの順である．
領域番号
$\sigma_x$
$\sigma_y$
$\tau_{xy}$
$\varepsilon_x$
$\varepsilon_y$
$\gamma_{xy}$
< ELEMENT STRESS & STRAIN >
ELEM.NO. SIG-X SIG-Y TAU-XY EPS-X EPS-Y GAM-XY
1 -0.18406E+03 0.42833E+03 0.20196E+03 -0.16727E-05 0.23078E-05 0.26255E-05
2 0.38602E+03 -0.20728E+04 -0.18426E+03 0.57984E-05 -0.10184E-04 -0.23954E-05
3 0.38602E+03 -0.20728E+04 0.18426E+03 0.57984E-05 -0.10184E-04 0.23954E-05
4 -0.18406E+03 0.42833E+03 -0.20196E+03 -0.16727E-05 0.23078E-05 -0.26255E-05
5 -0.65199E+03 -0.12644E+04 0.63409E+03 -0.50101E-06 -0.44816E-05 0.82432E-05
6 0.17275E+03 -0.38668E+04 0.27855E+03 0.83263E-05 -0.17931E-04 0.36212E-05
7 0.10753E+04 -0.17774E+04 0.16007E+03 0.83586E-05 -0.10184E-04 0.20809E-05
8 0.21775E+03 -0.43312E+04 0.41651E-05 0.94366E-05 -0.20132E-04 0.54147E-13
9 0.10753E+04 -0.17774E+04 -0.16007E+03 0.83586E-05 -0.10184E-04 -0.20809E-05
10 0.17275E+03 -0.38668E+04 -0.27855E+03 0.83263E-05 -0.17931E-04 -0.36212E-05
11 -0.65199E+03 -0.12644E+04 -0.63409E+03 -0.50101E-06 -0.44816E-05 -0.82432E-05
12 -0.15931E+04 -0.47475E+04 0.45770E+00 0.20088E-05 -0.18494E-04 0.59501E-08
13 -0.42182E+03 -0.48234E+04 0.69639E+02 0.74863E-05 -0.21124E-04 0.90531E-06
14 0.12123E+04 -0.41419E+04 0.14931E+03 0.13593E-04 -0.21210E-04 0.19411E-05
15 0.13368E+04 -0.38516E+04 -0.20428E-05 0.13593E-04 -0.20132E-04 -0.26556E-13
16 0.12123E+04 -0.41419E+04 -0.14931E+03 0.13593E-04 -0.21210E-04 -0.19411E-05
17 -0.42182E+03 -0.48234E+04 -0.69639E+02 0.74863E-05 -0.21124E-04 -0.90531E-06
18 -0.15931E+04 -0.47475E+04 -0.45770E+00 0.20088E-05 -0.18494E-04 -0.59501E-08
図7.18 領域内の応力とひずみの出力結果
一方，ペナルティ辺の数だけ降伏判定フラッグと単位面積当りの表面力が
図7.19 に示すように出力される．出力内容は左から以下のようになっている．降伏判定フラッグが0のとき弾性，1のとき降伏，2のとき除荷状態を表している．本例題では内部摩擦角を0としているため，破壊は垂直応力の如何に関わらず，せん断応力がせん断強度 $c$ に達した段階で生じている．すなわち，すべり破壊しているばねでは降伏判定フラッグがlとなり，せん断力の項が $0.3 \times 10^4$ と表示されている．
ペナルティ番号
降伏判定フラッグ
垂直応力
せん断応力
< STRESS >
SPRING.NO. YIELD SIG-N TAU-S
1 0 -0.95488E+02 -0.52586E+03
2 0 -0.52586E+03 0.95488E+02
3 0 -0.82512E+03 -0.64690E-05
4 0 0.43037E+03 0.20377E+03
5 0 -0.95488E+02 0.52586E+03
6 0 0.43037E+03 -0.20377E+03
7 0 -0.52586E+03 -0.95488E+02
8 0 0.95488E+02 0.52586E+03
9 0 -0.28679E+04 0.50405E+02
10 1 -0.23068E+04 0.30000E+04
11 0 -0.20168E+04 -0.16280E+04
12 0 -0.20168E+04 0.16280E+04
13 0 0.16010E+04 0.72414E-06
14 0 -0.28679E+04 -0.50405E+02
15 0 0.95488E+02 -0.52586E+03
16 1 -0.23068E+04 -0.30000E+04
17 1 -0.44104E+04 0.30000E+04
18 1 -0.33656E+04 -0.30000E+04
19 1 -0.65647E+02 0.30000E+04
20 1 -0.65647E+02 -0.30000E+04
21 1 -0.33656E+04 0.30000E+04
22 1 -0.44104E+04 -0.30000E+04
図7.19 ペナルティ辺の降伏判定フラッグと表面力の出力結果
次に，拘束辺の反力が
図7.20 のように出力される．ただし，出力内容は，下記のように左から，拘束番号，拘束する境界辺を持つ領域番号，拘束する境界辺の節点番号，拘束自由度方向，反力で，拘束自由度方向が1のとき，$x$ 方向，2のとき $y$ 方向を表している．
拘束番号
拘束領域番号
節点番号1
節点番号2
拘束方向
反力
< Reactopn:edge >
no. Elem.no. Node.no.1 Node.no.2 direc. Reaction
1 12 14 7 2 -0.92629E+02
2 12 11 14 1 0.17629E+02
3 14 12 11 1 -0.17629E+02
4 16 13 12 1 -0.17629E+02
5 18 15 13 1 0.17629E+02
6 18 10 15 2 0.92629E+02
図7.20 拘束辺の反力の出力結果
荷重制御によって，ペナルティ辺が降伏し，構造が崩壊機構条件を満たすと，最終的に連立方程式を解くことができなくなる．この時点で構造物は崩壊したとしてメッセージを出力して計算を終了させる．一方，計算誤差などの関係により，メカニズムに達していながら変位が求まってしまうこともある．しかし，それ以後の計算には意味が無いので，許容値以上の変位が生じていたならば
図7.21 に示すよう，計算された最大変位と許容値を出力し，計算を終了させる．ただし，局所破壊などのように，部分的にメカニズムを形成した場合も同様なメッセージが出力され計算が終了してしまうので注意が必要である．
********************************************
COLLAPS !! DMAX = 0.42280E+10 ( 0.10000E+04 )
********************************************
図7.21 変位が大きくなった場合のメカニズム・メッセージ
なお，本例題は，変位制御で解析しているため，
図7.21 のようなメッセージは出力されない．
以上の内容が，荷重ステップ分だけ繰り返して出力される．
