3DCAD3æœ‰é™ è¦ ç´ æ³..

長岡モノづくりアカデミー 3D-CAD/CAE コース
0）有限要素法の概略
’15.9
S
入力
形状、材料、荷重、拘束条件
形状
E,ν
要素剛性マトリックスK作成
Bマトリックス
Dマトリックス
全体剛性マトリックスK作成
f = K u
要素剛性マトリックスを対
応した行列に配置する
K の逆行列を求める。
変位　u = K -1 f を求める。
各要素の歪と応力を求め
る。
ε = B u, σ = D ε
必要に応じ、主応力、ミー
ゼス応力求める
E
1
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
f
１） 2次元形状を例に変位量
f により三角形123 が三角形1’2’3’ に変形
した。接点 1,2,3 の x 方向変位量は
u1  x'1  x1 

u 2  x'2  x 2 
u 3  x'3  x3 
Ω
(1.1)
さらに、 y 方向も同様に考えると要素におけ
る変形が表現できる。
今、変位量 u を便宜上 z 方向と考えると、
3点 (x1,y1,u1), (x2,y2,u2), (x3,y3,u3) で出来る
三角形は3次元空間での平面をなしてると仮
定すると、平面の式より以下のように書ける。
f
領域　Ω　の変形
1'
1
P'×
u  ax  by  c
(1.2)
この仮定は後で解るが、要素内の歪一定の
仮定と同じことである。
y
2
2'
P×
3'(x' ,y' )
3 3
3(x ,y )
3 3
変位量 (u3,v3)
x
2
長岡モノづくりアカデミー 3D-CAD/CAE コース
接点3点での座標と変位が求まっているとするなら
行列で書くと
u1  ax1  by1  c
u1   x1
u 2  ax 2  by 2  c u 2    x 2

u 3  ax3  by 3  c u 3   x3
’15.9
1'
1
P'×
y1 1 a 
 
y 2 1 b 
y 3 1 c 
y
(1.3)
2
2'
P×
3'(x' ,y' )
3 3
3(x ,y )
3 3
変位量 (u3,v3)
x
三角形要素の変形
座標と変位がわかっているので、平面を表す係数が求まる。
a   x1
  
b    x 2
c   x3
  
y1 1
y 2 1
y 3 1
1
u1 
 
u 2 
u 3 
 
a 
 y 2  y3
  1 
x3  x 2
b  

c  2  x 2 y 3  x3 y 2
 

この逆行列を実際に計算する。
2  x 2 y 3  x 2 y1  x1 y 3  x3 y 2  x3 y1  x1 y 2
で与えられ、Δは三角形の面積である。
y 3  y1
y1  y 2  u1 
 

x1  x3
x 2  x1  u 2 
x3 y1  x1 y 3 x1 y 2  x 2 y1 u 3 
(1.4)
3
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
これで係数が求まったので三角形要素内の任意の点Pでのｘ方向変位が求まる。
u  axP  byP  c
(1.5)
ｙ方向に関してもｘ方向と同じに考えると P点でのy方向変位が求まる。
v  a ' xP  b' yP  c'
(1.6)
ここで、係数は
a '
 y 2  y3
  1 
x3  x 2
b'  

c'  2  x 2 y 3  x3 y 2
 

y 3  y1
y1  y 2  v1 
 

x1  x3
x 2  x1  v 2 
x3 y1  x1 y 3 x1 y 2  x 2 y1 v3 
(1.4)'
とｘ方向のときと同様に求まる。
(1.5) (1.6) の２式より、三角形要素内の任意で座標を決めれば、変位が求
まる。
4
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
２）要素内部での変位量とひずみの関係
OABC の長方形が変位し、O’A’B’C’
に変形する。
x 方向の歪は
O' A'OA O' A' 'OA

OA
OA
u
u  dx  u
u
x


(1.7)
dx
x
x 
同様に y 方向は
O'C'  OC O' C ' 'OC

OC
OC
v
v  dy  v
v
y


(1.8)
dy
y
y 
せん断歪は
xy  tan   tan  

v u

x y
A' A' ' C ' C ' '

O' A' ' O' C ' '
(1.9)
5
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
ここで先に求めた、下の要素内座標と変位の式を用いると
u  axP  byP  c
v  a' xP  b' yP  c'
(1.5)
(1.6)
各歪が変位より求めることができ、その値は要素内では一定になっている。
さらに、係数と変位の関係式（１．４）式を合わせ考えると
 u




x
 a 
x  
 
   v


y

  
  b' 
xy   y
 

   v u  a 'b 
  
 x y 
 y 2  y3
1 

0
2 
 x3  x 2
0
x3  x 2
y 2  y3
(1.10)
y 3  y1
0
x1  x3
0
x1  x3
y 3  y1
y1  y 2
0
x 2  x1
u1 
v1 
0  
u 2 

x 2  x1   
v2
y1  y 2   
u 3 
 
v3 
(1.11)
6
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
これらの式が意味するところは
与えられた領域を小さい要素に分割し、
その要素の各３点の座標を決めてやり、各要素内での歪を
一定と考えると（三角形要素を使うことを意味する）
各要素内の変位は要素内の座標を与えると
座標点がわかると作ることができる行列
（これをBマトリックスと呼ぶ）と
各要素の接点での変位ベクトルの掛け算で
要素内の歪を求めることができる。つまり変位が求められれ
ば歪、ひいては応力が求められる。
7
長岡モノづくりアカデミー 3D-CAD/CAE コース
 y 2  y3
1 
B
0
2 
 x3  x 2
0
x3  x 2
y 2  y3
y 3  y1
0
x1  x3
0
x1  x3
y 3  y1
y1  y 2
0
x 2  x1
0 
x 2  x1 
y1  y 2 
’15.9
(1.12)
として、要素内の変位 u とひずみ ε 、応力 σ とひずみ ε の関係を書く
と、次のように書ける。
ε  Bu


σ  Dε  DBu 
(1.13)
ここで Dは平面応力状態なので
1 
 x
0
x 
 


E
 
 
 1
σ  y  
0  y   Dε
xy  (1  )(1  ) 0 0 1   xy 
 
 

2 
こう書ける理由は別に述べる。
8
長岡モノづくりアカデミー 3D-CAD/CAE コース
F(σ)
  DBu
F
A
dU
U 
u
du
U
O
’15.9
u
F
dU
du
  Bu
A
W  Fu
W
u(ε)
歪を要素全体で積分した
歪エネルギー
O
δu
u
u
外力によりなされる仕事
今弾性域（線形域内）で考えている。Ａ点で釣り合っている。今荷重点ＡにＦが
作用し釣り合いを保っている。そのとき要素内は歪によりＵの歪エネルギーが
たまっている。今この接点に作用する力Ｆで仮想変位δu 変位した。そのとき、外
力がする仕事Ｗはすべて歪エネルギーとして要素内部に蓄えられる。
 dU

U  W  
 F u  0
 du

よって
dU
F
du
9
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
ここで、要素全体の歪エネルギーは弾性範囲内を考えているので
1 T
1
1 T T T
T
U   σ εd   ( DBu) ( Bu)d   u B D Bud
2
2
2


1 T
u Kud

2
ここで、
K  BT DT B
さらに、Ｋは対称行列、よってＵは2次形式であるので
dU  U U U U U U 

,
,
,
,
,
   Kud
du  u1 v1 u 2 v 2 u 3 v3 

T
よって
 Kud  F
と成る。

10
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
さてここで、
 Kud
は要素の座標が決まればＫは定数である。

 Kud  Ku  d  Kut


ここで第2項の積分は三角形要素全体での積分でありつまり要素の体積
である。よって Δ は三角形の面積、ｔは三角形要素の厚さで表される。
よって
Kt  K
と再度置き直す。すると各要素において
Ku  F
u  K 1F
u  ( BT DT B) 1 F
と外部からの接点外力に応じた接点変位が
求まる。
11
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
２）三角形要素を用い、片持ち梁の問題を解く
さて、ここまでの話は各要素内での一般的な話であった。
以降は、片持ち梁を例に
Ａ）
Ｂ)
Ｃ）
Ｄ）
Ｅ）
Ｆ）
要素剛性マトリックスを作り
全体の剛性マトリックスにまとめ
外力より、拘束点以外の接点変位ベクトルを求め
全変位ベクトルから固定点の節点反力を求める。
要素における変位ベクトルから各要素内の歪を求める。
歪から要素内の応力を求める。
さらにここまで求まると
Ｇ）要素内の応力より主応力、ミーゼス応力等を求めることが出来る。
それにより、材料の安全性の検証可能となる。
この件に関しては別のスライドで説明する。
12
長岡モノづくりアカデミー 3D-CAD/CAE コース
f=(0,1)
y
1
5
6
図のような、一辺が１の三角形要素4
個で出来た片持ち梁を考え、節点4に
外力（１，０）が作用。
節点1，6が固定されている。
4
③
④
要素①
節点1
0
’15.9
②
2
1
3
2
本来、右の式ではあるが
今回は計算を楽にするため
に
 1 0.5 0
D  0.5 1 0
 0
0 1
x


1

0
x 
x 


E
 
 
 1
σ  y  
0  y   Dε
xy  (1  )(1  ) 0 0 1   xy 
 

 
2 

(2.1)
とおいて、計算を進める。
13
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
有限要素法を用いて片持ち梁を計算するとき、未知数は変位である。
歪、応力、節点反力等は、節点変位から計算される。各要素において
Ku  F
(2.2)
u  K 1F
(2.2)'
u  ( BT DT B ) 1 F
(2.3)
で与えられることは既に述べた。
 y 2  y3
1 
B
0

2
 x3  x 2
0
x3  x 2
y 2  y3
y 3  y1
0
x1  x3
0
x1  x3
y 3  y1
y1  y 2
0
x 2  x1
0 
x 2  x1 
y1  y 2 
(1.12)
であり、三角形要素の面積は (1/2)(底辺)(高さ)＝1/2 であるから、
14
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
要素①のBマトリックスは
0
y 5  y1
 y 2  y5
1 
B① 
0
x5  x 2
0
1
2   x5  x 2 y 2  y 5 x1  x5
2
  1 0 1 0 0 0
  0 0 0  1 0 1
 0  1  1 1 1 0
0
x1  x5
y 5  y1
y1  y 2
0
x 2  x1
0 
x 2  x1 
y1  y 2 
(2.11)
要素②のBマトリックスは同様に、座標の差で与えられているのでまったく同じに
なる。
  1 0 1 0 0 0
B ①  B ②   0 0 0  1 0 1
 0  1  1 1 1 0
(2.12)
15
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
要素③のBマトリックスは
0
y5  y 2
0
 y 4  y5
1 
B③ 
0
x5  x 4
0
x 2  x5
1
2   x5  x 4 y 4  y 5 x 2  x5 y 5  y 2
2
 0 0 1 0 1 0 
  0  1 0 0 0 1 
(2.13)
 1 0 0 1 1  1
y2  y4
0
x4  x2
0 
x 4  x 2 
y 2  y 4 
要素④のBマトリックスは同様に、座標の差で与えられているのでまったく同じに
なり、
 0 0 1 0 1 0 
B ③  B ④   0  1 0 0 0 1 
 1 0 0 1 1  1
(2.14)
16
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
各要素におけるＢマトリックスとＤマトリックスが求められたので、
各要素ごとの要素剛性マトリックスＫが求められる。
F  Ku
のＫは要素①、②では
 1 0 0 
 0 0  1
  1 0.5 0  1 0 1 0 0 0

 1 0  1 


K ①  K ②  BT DT B  
 0.5 1 0  0 0 0  1 0 1
 0 1 1  0
0 1  0  1  1 1 1 0
 0 0 1 


 0 1 0 
0.5 0  0.5
1
0
 1

 0
0
1

1

1
1


 1
 1.5  1 0.5 
2
1


1

1
2
5
.
1

1

5
.
0


 0
0 
1
1
1 1


1
0
1

5
.
0
0
5
.
0



(2.15)
17
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
同じく要素③、④では
 0 0  1
 0 1 0 
  1 0.5 0  0 0 1 0  1 0 


1 0 0 


K ③  K ④  BT DT B  
 0.5 1 0  0  1 0 0 0 1 
 0 0 1  0
0 1  1 0 0 1 1  1
 1 0 1  


 0 1  1
1 
1 1
0
0
1
0
 1 
0.5
 05 0
1

 0  0.5
0.5 
1
0
(2.16)


5
.
0

1
1
0
0
1



 1 0.5
 1.5 
2
1 1


2 
0.5  1  1.5
1
 1
以上は各要素ごとの話であったが、これを全体の剛性マトリックスにまとめる。
つまり、Ｆ＝Ｋｕのｕ、Ｆは節点1から6まで全ての項が並ぶ。
よって、Ｋは12×12 の全体の剛性マトリックスとなる。このキングサイズマト
リックスに各要素の項を入れ込むのである。
18
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
例えば、要素①の場合
0
1
0.5
 f 1
 1
 g1 
 0
1
1
1
 

 f 2
 1
1
2
 1.5
 

g
2
2
 
 0.5  1  1.5
 f 3
 0
0
0
0
 

0
0
0
 g 3   F  Ku   0
 f 4
 0
0
0
0
 

0
0
0
 g 4
 0
 f 5
 0
1 1
1
 

 g5
 0.5 0
0.5
1
 

0
0
0
 f 6
 0
 g 6 
 0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0 0  0.5
0 1
0
0  1 0.5
0 1
1
0 0
0
0 0
0
0 0
0
0 0
0
0 1
0
0 0
1
0 0
0
0 0
0
0
0
0
0
0
0
0
0
0
0
0
0
0   u1 
0  v1 
0 u 2 
 
0  v 2 
0  u 3 
 
0  v 3 
0 u 4 
 
0  v 4 
0  u 5 
 
0  v 5 
 
0 u 6 
0  v 6 
先程、要素剛性マトリックスのところで求めたＫの要素を要素①を
作っている節点１，２，５が関係するところに配置しなおす。
19
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
もうひとつ例を示すと、要素③の場合
要素剛性マトリックスのところで求めたＫの要素を要素③を作ってい
る節点２，４，５が関係するところに配置しなおす。
 f 1
0
 g1 
0
 

 f 2
0
 

g
2
 
0
 f 3
0
 

g
3
   F  Ku  0
 f 4
0
 

 g 4
0
 f 5
0
 

 g5
0
 

 f 6
0
 g 6 
0
0 0
0
0 0
0
0 1
0
0 0
1
0 0
0
0 0
0
0 0  0.5
0 1
0
0  1 0.5
0 1
1
0 0
0
0 0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
1 1
1
0  0.5 0
0.5
1
0
0
0
0
0
0
0
0
0
0
0
1
0
1
0.5
0
0
1
1
1
0 1
1
2
 1.5
0 0.5  1  1.5
2
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0   u1 
0  v1 
0 u 2 
 
0  v 2 
0  u 3 
 
0  v 3 
0 u 4 
 
0  v 4 
0  u 5 
 
0  v 5 
 
0 u 6 
0  v 6 
20
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
このように、①、②、③、④に対する要素剛性マトリックスを全て配置して全体の
剛性マトリックスを作る。
K  K①  K②  K③  K④
求
め
る
反
力
よって、固定点や外部からの荷重等を考慮すると
0
1
0.5
0
0
0
0
0
 1.5  1
1  0 
 f 1  2
 g1   0
2
1
1
0
0
0
0
 1.5
0
0.5
 1   0 
  
 0   1
1
4
 1.5  1
0.5
0
 1.5  2 1.5
0
0  u 2 
  
 
0
0
.
5

1

1
.
5
4
1

1

1
.
5
0
1
.
5

2
0
0
  
 v2 
0  0
0
1
1
2
 1.5  1
0.5
0
0
0
0  u 3 
  
 
0
0.5
 1  1.5
2
1
1
0
0
0
0   v3 
0 0
0  0
0
0
 1.5  1
1
2
0
1
0.5
0
0  u 4 
  
 
1
0
0

1
.
5
0
0
.
5

1
0
2
1

1
0
0
  
 v4 
0  0
 1.5  2 1.5
0
0
1
1
4
 1.5  1
0.5  u 5 
  
 
 0   1.5
0
1.5
2
0
0
0.5
 1  1.5
4
1
 1   v5 
  
 
0.5
0
0
0
0
0
0
1
1
2
 1.5  0 
 f 6   1
 g 6   1
1
0
0
0
0
0
0
0.5
 1  1.5
2   0 
荷重
固
定
条
件
この部分をＫｐとし、逆行列もとめる
と、既知の力から変位が求められる。
21
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
つまり、以下のように表せる。
 1.5  1
0.5
0
 1.5  2 1.5  u 2 
0 
 4
0 
 1.5
4
1
 1  1.5
0
1.5
 2   v 2 
 

0 
 1
1
2
 1.5  1
0.5
0
0  u 3 
 

 
 1  1.5
2
1
1
0
0   v3 
0  Kpu   0.5
0 
 0
 1.5  1
1
2
0
1
0.5  u 4 
 

 
0
0.5
1
0
2
1
 1  v4 
1
 1.5
0 
  2 1.5
0
0
1
1
4
 1.5 u 5 
 

 
2
0
0
0.5
 1  1.5
4   v5 
0
 1.5
u 2 
0.99 0.59
v2 
0.59 1.22
 

u 3 
0.96 0.55
 

v
3

1
   Kp F  1.07 1.69
u 4 
0.44 0.01
 

v
4
 
1.04 1.65
u 5 
0.20  0.32
 

v
5
 
0.20 0.68
0.96 1.07
0.44
1.04
0.55 1.69
0.01
1.65
2.17 1.95
0.57
1.15
1.96 4.85  0.89 3.74
0.57  0.89 1.67  0.76
1.15 3.74  0.76 3.86
0.21  0.81 0.71  0.81
0.21 1.19  0.29 1.19
逆行列を求めるプログラ
ムは在るので
0.20
0.20  0  1.04 
 0.32 0.68  0  1.65 
0.21
0.21  0  1.15 
  

 0.81 1.19  0  3.74 

0.71  0.29 0  0.76
  

 0.81 1.19  1  3.86 
0.76  0.24 0   0.81
  

 0.24 0.76  0  1.19 
(2.28)
と、変位が求まる。
22
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
以上より、変位が計算されたので、
F  Ku
0
1
0.5
0
0
0
0
0
 1.5  1
1  0 
 2   2
 0.60   0
2
1
1
0
0
0
0
 1.5
0
0.5
 1   0 

 
 0   1
1
4
 1.5  1
0.5
0
 1.5  2 1.5
0
0   1.04 

 


0
0
.
5

1

1
.
5
4
1

1

1
.
5
0
1
.
5

2
0
0
1
.
65

 


 0   0


0
1
1
2
 1.5  1
0.5
0
0
0
0
1.15 

 


0
0.5
 1  1.5
2
1
1
0
0
0
0   3.74 
 0  0
 0   0
0
0
 1.5  1
1
2
0
1
0.5
0
0   0.76

 


0
 1.5
0
0.5
1
0
2
1
1
0
0   3.86 
 1   0
 0   0
 1.5  2 1.5
0
0
1
1
4
 1.5  1
0.5    0.81

 


 0   1.5
0
1.5
2
0
0
0.5
 1  1.5
4
1
 1   1.19 

 


2

1
0
.
5
0
0
0
0
0
0

1
1
2

1
.
5
0

 


 1.6  1


1
0
0
0
0
0
0
0.5
 1  1.5
2   0 
(2.29)
以上のように、変位から節点１，６における固定点反力がｘ、ｙそれぞれの方向別
に求めることが出来る。
23
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
このように、各接点での反力と、各要素接点での変位が求められたので
ｘ、ｙ方向での、要素内の歪、要素内での応力が求められる。
ε  Bu


σ  Dε  DBu 
(1.13)
例えば、要素①での歪と応力は
 u1 
 0 
 v1 
 0 
    1 0 1 0 0 0 
  1.04 
 x 
u 2 
 1.04  

ε1   y   B①     0 0 0  1 0 1 
   0.46
 v 2   0  1  1 1 1 0  1.65    0.2 
xy 
  0.81 

u 5  
 


v
5
1
.
19


 

0.5 0  1.04   0.81 
x 
 1


σ1  y   Dε1  0.5 1 0  0.46   0.06 


 0
   0.2 
xy 
0
1
  0.2


(2.30)
(2.31)
以上のように、要素内での歪と応力がｘ、ｙ方向に対して求めることが出来る。
この応力を使い主応力やミーゼス応力を計算し、材料が応力に耐えうるかどうか
の検証が可能となる。本件は別に検討する。
24
長岡モノづくりアカデミー 3D-CAD/CAE コース
’15.9
25

Download Report