School of Civil, Environmental and Mining Engineering
Finite Element Theory and Practice
2D CST element: examples
Example 1
Example 1: Evaluate the strain and stress components at any point within a 3 node triangular plane stress
element. The coordinates of the 3 nodes & corresponding nodal displacements are given below. Take E =
100GPa & Poisson’s ratio: ν = 0.25
𝑥1 = 1.2𝑚, 𝑦1 = 0.5𝑚, 𝑥2 = 1.8𝑚, 𝑦2 = 0.8𝑚, 𝑥3 = 1.6𝑚, 𝑦3 = 1.2𝑚
𝑢1 = 5.6𝑚𝑚, 𝑣1 = 4.3𝑚𝑚, 𝑢2 = 6.4𝑚𝑚, 𝑣2 = 4.6𝑚𝑚, 𝑢3 = 6.1𝑚𝑚, 𝑣3 = 5.1𝑚𝑚
Strain vector for 2D plane stress/strain:
𝜀𝑥 𝜕𝑢/𝜕𝑥
𝜀 = 𝜀𝑦 = 𝜕𝑣/𝜕𝑦
𝛾𝑥𝑦 𝜕𝑢/𝜕𝑦 + 𝜕𝑣/𝜕𝑥
where displacement fields 𝑢 and 𝑣 are approximated as
𝑢 = 𝑁1 𝑢1 + 𝑁2 𝑢2 + 𝑁3 𝑢3
𝑣 = 𝑁1 𝑣1 + 𝑁2 𝑣2 + 𝑁3 𝑣3
First component of strain vector, normal strain 𝜀𝑥 , is:
𝜕𝑢 𝜕𝑁1 𝜕𝑁2 𝜕𝑁3
𝜀𝑥 = = 𝑢 + 𝑢 + 𝑢
𝜕𝑥 𝜕𝑥 1 𝜕𝑥 2 𝜕𝑥 3
Example 1
First term, normal strain 𝜀𝑥 , can be approximated as:
𝜕𝑢 𝜕𝑁1 𝜕𝑁2 𝜕𝑁3
𝜀𝑥 = = 𝑢 + 𝑢 + 𝑢
𝜕𝑥 𝜕𝑥 1 𝜕𝑥 2 𝜕𝑥 3
Shape functions and their derivatives
1 1 1
𝑁1 = 𝑎 + 𝑏1 𝑥 + 𝑐1 𝑦 𝑁2 = 𝑎 + 𝑏2 𝑥 + 𝑐2 𝑦 𝑁3 = 𝑎 + 𝑏3 𝑥 + 𝑐3 𝑦
2Δ 1 2Δ 2 2Δ 3
𝜕𝑁1 𝑏1 𝜕𝑁2 𝑏2 𝜕𝑁3 𝑏3
= = =
𝜕𝑥 2Δ 𝜕𝑥 2Δ 𝜕𝑥 2Δ
where
𝑎1 = 𝑥2 𝑦3 − 𝑦2 𝑥3 = 1.8 × 1.2 − 0.8 × 1.6 = 0.88𝑚2 𝑏1 = 𝑦2 − 𝑦3 = 0.8 − 1.2 = −0.4𝑚 𝑐1 = 𝑥3 − 𝑥2 = 1.6 − 1.8 = −0.2𝑚
𝑎2 = 𝑥3 𝑦1 − 𝑦3 𝑥1 = 1.6 × 0.5 − 1.2 × 1.2 = −0.64𝑚2 𝑏2 = 𝑦3 − 𝑦1 = 1.2 − 0.5 = 0.7𝑚 𝑐2 = 𝑥1 − 𝑥3 = 1.2 − 1.6 = −0.4𝑚
𝑎3 = 𝑥1 𝑦2 − 𝑦1 𝑥2 = 1.2 × 0.8 − 0.5 × 1.8 = 0.06𝑚2 𝑏3 = 𝑦1 − 𝑦2 = 0.5 − 0.8 = −0.3𝑚 𝑐3 = 𝑥2 − 𝑥1 = 1.8 − 1.2 = 0.6𝑚
2Δ = 𝑎1 + 𝑎2 + 𝑎3 = 0.88 − 0.64 + 0.06 = 0.3𝑚2
Normal strain 𝜀𝑥 :
𝑏1 𝑏2 𝑏3 −0.4 × 5.6 × 10−3 + 0.7 × 6.4 × 10−3 + −0.3 × 6.1 × 10−3
𝜀𝑥 = 𝑢1 + 𝑢2 + 𝑢 =
2Δ 2Δ 2Δ 3 0.3
= 1.366 × 10−3
Example 1
1 𝜕𝑁1 𝑏1
𝑁1 = 𝑎 + 𝑏1 𝑥 + 𝑐1 𝑦 =
Second term, normal strain 𝜀𝑦 , can be obtained in a similar way: 2Δ 1 𝜕𝑥 2Δ
𝜕𝑣 𝜕𝑁1 𝜕𝑁2 𝜕𝑁3 𝑐1 𝑐2 𝑐3 1 𝜕𝑁2 𝑏2
𝜀𝑦 = = 𝑣 + 𝑣 + 𝑣 = 𝑣1 + 𝑣2 + 𝑣 𝑁2 = 𝑎 + 𝑏2 𝑥 + 𝑐2 𝑦 =
𝜕𝑦 𝜕𝑦 1 𝜕𝑦 2 𝜕𝑦 3 2Δ 2Δ 2Δ 3 2Δ 2 𝜕𝑥 2Δ
1 𝜕𝑁3 𝑏3
−0.2 × 4.3 × 10−3 + −0.4 × 4.6 × 10−3 + 0.6 × 5.1 × 10−3 𝑁3 = 𝑎 + 𝑏3 𝑥 + 𝑐3 𝑦 𝜕𝑥
=
2Δ
= 2Δ 3
0.3
= 1.2 × 10−3 𝑎1 = 𝑥2 𝑦3 − 𝑦2 𝑥3 = 1.8 × 1.2 − 0.8 × 1.6 = 0.88𝑚2
𝑎2 = 𝑥3 𝑦1 − 𝑦3 𝑥1 = 1.6 × 0.5 − 1.2 × 1.2 = −0.64𝑚2
Shear strain 𝛾𝑥𝑦 :
𝑎3 = 𝑥1 𝑦2 − 𝑦1 𝑥2 = 1.2 × 0.8 − 0.5 × 1.8 = 0.06𝑚2
𝜕𝑢 𝜕𝑣 𝜕𝑁1 𝜕𝑁2 𝜕𝑁3 𝜕𝑁1 𝜕𝑁2 𝜕𝑁3
𝛾𝑥𝑦 = + = 𝑢1 + 𝑢2 + 𝑢3 + 𝑣1 + 𝑣2 + 𝑣 2Δ = 𝑎1 + 𝑎2 + 𝑎3 = 0.88 − 0.64 + 0.06 = 0.3𝑚2
𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜕𝑦 𝜕𝑦 𝜕𝑥 𝜕𝑥 𝜕𝑥 3
𝑏1 = 𝑦2 − 𝑦3 = 0.8 − 1.2 = −0.4𝑚
𝑐1 𝑐2 𝑐3 𝑏1 𝑏2 𝑏3
= 𝑢 + 𝑢 + 𝑢 + 𝑣 + 𝑣 + 𝑣
2Δ 1 2Δ 2 2Δ 3 2Δ 1 2Δ 2 2Δ 3 𝑏2 = 𝑦3 − 𝑦1 = 1.2 − 0.5 = 0.7𝑚
𝑏3 = 𝑦1 − 𝑦2 = 0.5 − 0.8 = −0.3𝑚
−0.2 × 5.6 × 10−3 + −0.4 × 6.4 × 10−3 + 0.6 × 6.1 × 10−3
+ −0.4 × 4.3 × 10−3 + 0.7 × 4.6 × 10−3 + −0.3 × 5.1 × 10−3 𝑐1 = 𝑥3 − 𝑥2 = 1.6 − 1.8 = −0.2𝑚
=
0.3 𝑐2 = 𝑥1 − 𝑥3 = 1.2 − 1.6 = −0.4𝑚
= −0.166 × 10−3 𝑐3 = 𝑥2 − 𝑥1 = 1.8 − 1.2 = 0.6𝑚
Example 1
𝜀𝑥 1.366
Strain vector obtained:
𝜀 = 𝜀𝑦 = 1.2 10−3
𝛾𝑥𝑦 −0.166
Stress-strain relationship: 𝜎 = 𝐷 𝜀
𝐸 1 𝜈 0 1 0.25 0 1 0.25 0
100
𝐷 = 𝜈 1 0 = 0.25 1 0 = 106.667 0.25 1 0 GPa
1 − 𝜈 2 0 0 (1 − 𝜈)/2 1 − 0.252 0 0 (1 − 0.25)/2 0 0 0.375
(𝐸 = 100GPa, 𝜈 = 0.25)
Stress vector:
x 1 0.25 0 1.366 177.707
−3
y = 106.667 0.25 1 0 1.2 10 GPa = 164.427 MPa
0 0.375 −0.166 −6.640
xy 0
Example 1
1 𝜕𝑁1 𝑏1
Alternatively the strain-displacement relationship can be used: 𝑁1 = 𝑎 + 𝑏1 𝑥 + 𝑐1 𝑦
2Δ 1 =
𝜕𝑥 2Δ
1 𝑏1 0 𝑏2 0 𝑏3 0 𝜕𝑁2 𝑏2
1
𝜀 = 𝐵 𝑋 𝐵 = 0 𝑐1 0 𝑐2 0 𝑐3 𝑁2 = 𝑎 + 𝑏2 𝑥 + 𝑐2 𝑦 =
2Δ 2Δ 2 𝜕𝑥 2Δ
𝑐1 𝑏1 𝑐2 𝑏2 𝑐3 𝑏3
1 𝜕𝑁3 𝑏3
𝑁3 = 𝑎 + 𝑏3 𝑥 + 𝑐3 𝑦 =
5.6 2Δ 3 𝜕𝑥 2Δ
4.3
6.4
𝑋 = 103 × 𝑎1 = 𝑥2 𝑦3 − 𝑦2 𝑥3 = 1.8 × 1.2 − 0.8 × 1.6 = 0.88𝑚2
4.6
6.1 𝑎2 = 𝑥3 𝑦1 − 𝑦3 𝑥1 = 1.6 × 0.5 − 1.2 × 1.2 = −0.64𝑚2
5.1 𝑎3 = 𝑥1 𝑦2 − 𝑦1 𝑥2 = 1.2 × 0.8 − 0.5 × 1.8 = 0.06𝑚2
5.6 × 10−3 2Δ = 𝑎1 + 𝑎2 + 𝑎3 = 0.88 − 0.64 + 0.06 = 0.3𝑚2
4.3 × 10−3
1 −0.4 0 0.7 0 −0.3 0
6.4 × 10−3 𝑏1 = 𝑦2 − 𝑦3 = 0.8 − 1.2 = −0.4𝑚
𝜀 = 0 −0.2 0 −0.4 0 0.6
0.3 4.6 × 10−3 𝑏2 = 𝑦3 − 𝑦1 = 1.2 − 0.5 = 0.7𝑚
−0.2 −0.4 −0.4 0.7 0.6 −0.3
6.1 × 10−3
5.1 × 10−3 𝑏3 = 𝑦1 − 𝑦2 = 0.5 − 0.8 = −0.3𝑚
1.3667 𝑐1 = 𝑥3 − 𝑥2 = 1.6 − 1.8 = −0.2𝑚
= 103 × 1.2 𝑐2 = 𝑥1 − 𝑥3 = 1.2 − 1.6 = −0.4𝑚
−0.16667 𝑐3 = 𝑥2 − 𝑥1 = 1.8 − 1.2 = 0.6𝑚
Example 2
Example 2: The 3 node plane stress triangular element as shown below (a = 1.0m, b = 0.5m, c = 2m) is subjected
uniformly distributed body force of 200N/m3 acting along the direction of y. Evaluate the nodal load vector of the
element taking unit thickness (t =1.0m).
Nodal loads due to distributed body force:
𝑁1 0 𝑁1 𝑞𝑥
0 𝑁1 𝑁1 𝑞𝑦
𝑁 0 𝑞𝑥 𝑁2 𝑞𝑥
𝑃 = 𝑡න 𝑁 𝑇 𝑞 𝑑𝐴 = 𝑡 න 2 𝑑𝐴 = 𝑡 න 𝑁2 𝑞𝑦 𝑑𝐴
0 𝑁2 𝑞𝑦
𝑁3 0 𝑁3 𝑞𝑥
0 𝑁3 𝑁3 𝑞𝑦
where shape functions are:
1 1 1
𝑁1 = 𝑎 + 𝑏1 𝑥 + 𝑐1 𝑦 𝑁2 = 𝑎 + 𝑏2 𝑥 + 𝑐2 𝑦 𝑁3 = 𝑎 + 𝑏3 𝑥 + 𝑐3 𝑦
2Δ 1 2Δ 2 2Δ 3
Area integration is needed; before that, we need all coefficients 𝑎𝑖 , 𝑏𝑖 , 𝑐𝑖 of the shape functions.
Example 2
Nodal loads due to distributed body force:
𝑁1 𝑞𝑥 1
𝑁1 = 𝑎 + 𝑏1 𝑥 + 𝑐1 𝑦
𝑁1 𝑞𝑦 2Δ 1
𝑁2 𝑞𝑥 1
𝑃 = 𝑡 න 𝑁 𝑞 𝑑𝐴 𝑁2 = 𝑎 + 𝑏2 𝑥 + 𝑐2 𝑦
2 𝑦 2Δ 2
𝑁3 𝑞𝑥 1
𝑁3 𝑞𝑦 𝑁3 = 𝑎 + 𝑏3 𝑥 + 𝑐3 𝑦
2Δ 3
We need all coefficients 𝑎𝑖 , 𝑏𝑖 , 𝑐𝑖 of the shape functions. In local coordinate system:
𝑥1 = −𝑏 = −0.5𝑚 𝑥2 = 𝑎 = 1𝑚, 𝑥3 = 0,
𝑦1 = 0, 𝑦2 = 0, 𝑦3 = 𝑐 = 2𝑚
The coefficients:
𝑎1 = 𝑥2 𝑦3 − 𝑦2 𝑥3 = 1 × 2 − 0 × 0 = 2𝑚2 𝑏1 = 𝑦2 − 𝑦3 = 0 − 2 = −2𝑚 𝑐1 = 𝑥3 − 𝑥2 = −1𝑚
𝑎2 = 𝑥3 𝑦1 − 𝑦3 𝑥1 = 0 × 0 − 2 × (−0.5) = 1𝑚2 𝑏2 = 𝑦3 − 𝑦1 = 2 − 0 = 2𝑚 𝑐2 = 𝑥1 − 𝑥3 = −0.5𝑚
𝑎3 = 𝑥1 𝑦2 − 𝑦1 𝑥2 = (−0.5) × 0 − 0 × 1 = 0 𝑐3 = 𝑥2 − 𝑥1 = 1.5𝑚
𝑏3 = 𝑦1 − 𝑦2 = 0 − 0 = 0
2Δ = 𝑎1 + 𝑎2 + 𝑎3 = 2 + 1 + 0 = 3𝑚2
Example 2
Shape functions: 𝑎1 = 2𝑚2 𝑎2 = 1𝑚2 𝑎3 = 0
1 2 − 2𝑥 − 𝑦
𝑁1 = 𝑎 + 𝑏1 𝑥 + 𝑐1 𝑦 = 𝑏1 = −2𝑚 𝑏2 = 2𝑚 𝑏3 = 0
2Δ 1 3
𝑐1 = −1𝑚 𝑐2 = −0.5𝑚 𝑐3 = 1.5𝑚
1 1 + 2𝑥 − 0.5𝑦
𝑁2 = 𝑎 + 𝑏2 𝑥 + 𝑐2 𝑦 = 2Δ = 3𝑚2
2Δ 2 3
1 1.5𝑦 𝑞𝑥 = 0, 𝑞𝑦 = 200𝑁/𝑚3 𝑡 = 1𝑚
𝑁3 = 𝑎 + 𝑏3 𝑥 + 𝑐3 𝑦 =
2Δ 3 3
Load vector:
0
𝑁1 𝑞𝑥 0 2 න𝑑𝐴 − 2 න𝑥𝑑𝐴 − න𝑦𝑑𝐴
𝑁1 𝑞𝑦 2 − 2𝑥 − 𝑦
0
𝑁2 𝑞𝑥 1 × 200 0 200
𝑃 = 𝑡න 𝑁 𝑞 = න 𝑑𝐴 =
2 𝑦 3 1 + 2 − 0.5𝑦 3 න𝑑𝐴 + 2 න𝑥𝑑𝐴 − 0.5 න𝑦𝑑𝐴
𝑁3 𝑞𝑥 0
0
𝑁3 𝑞𝑦 1.5𝑦
1.5 න𝑦𝑑𝐴
Looking at the integrations in the above load vector, we need to evaluate the following terms:
න𝑑𝐴 = Δ = 3/2 = 1.5𝑚2 න𝑥𝑑𝐴 = ? න𝑦𝑑𝐴 = ?
Example 2
𝑎 = 1𝑚
Load vector:
0 𝑏 = 0.5𝑚
2 න𝑑𝐴 − 2 න𝑥𝑑𝐴 − න𝑦𝑑𝐴
𝑐 = 2𝑚
0
200
𝑃 =
3 න𝑑𝐴 + 2 න𝑥𝑑𝐴 − 0.5 න𝑦𝑑𝐴
0
1.5 න𝑦𝑑𝐴 𝑚! 𝑛!
න𝑥 𝑚 𝑦 𝑛 𝑑𝐴 = 𝑐 𝑛+1 𝑎𝑚+1 − −𝑏 𝑚+1
(𝑚 + 𝑛 + 2)!
Evaluation of terms, using the provided formula:
1 0 0+1 1+1 1+1
1! 0! 𝑐 𝑎2 − 𝑏 2 2× 1 2
− 0.5 2
න𝑥𝑑𝐴 = න𝑥 𝑦 𝑑𝐴 = 𝑐 𝑎 − −𝑏 = = = 0.25𝑚3
(0 + 1 + 2)! 6 6
0 1 1+1 0+1 0+1
0! 1! 𝑐2 𝑎 + 𝑏 2 2
× 1 + 0.5
න𝑦𝑑𝐴 = න𝑥 𝑦 𝑑𝐴 = 𝑐 𝑎 − −𝑏 = = = 1.0𝑚3
(1 + 0 + 2)! 6 6
Load vector:
0 0 0
2 × 1.5 − 2 × 0.25 − 1.0 1.5 100
200 0 200 0 0
𝑃 = = = 𝑘𝑁
3 0.15 + 2 × 0.25 − 0.5 × 1.0 3 1.5 100
0 0 0
1.5 × 1.0 1.5 100
