Transcript Lecture 6
MECH593 Introduction to Finite
Element Methods
Finite Element Analysis of 2-D Problems
Dr. Wenjing Ye
2-D Discretization
Common 2-D elements:
2-D Model Problem with Scalar Function
- Heat Conduction
• Governing Equation
T ( x, y) T ( x, y)
Q( x, y) 0
x
x y
y
• Boundary Conditions
Dirichlet BC:
Natural BC:
Mixed BC:
in W
Weak Formulation of 2-D Model Problem
• Weighted - Integral of 2-D Problem ---- T ( x, y ) T ( x, y )
W w x x y y Q( x, y )dA 0
• Weak Form from Integration-by-Parts ---- w T w T
0
wQ( x, y ) dxdy
x x y y
W
T
T
w
dxdy
dxdy w
x
x
y
y
W
W
Weak Formulation of 2-D Model Problem
• Green-Gauss Theorem ----Ω
Ω
𝜕
𝜕𝑇
𝜅𝑤
𝑑𝑥𝑑𝑦 =
𝜕𝑥
𝜕𝑥
𝜕
𝜕𝑇
𝜅𝑤
𝑑𝑥𝑑𝑦 =
𝜕𝑦
𝜕𝑦
Γ
Γ
𝜕𝑇
𝜅𝑤
𝑛 𝑑𝑠
𝜕𝑥 𝑥
𝜕𝑇
𝜅𝑤
𝑛𝑦 𝑑𝑠
𝜕𝑦
where nx and ny are the components of a unit vector,
which is normal to the boundary G.
n nx i n y j i cos j sin
Weak Formulation of 2-D Model Problem
• Weak Form of 2-D Model Problem ----𝜕𝑤
𝜕𝑇
𝜕𝑤
𝜕𝑇
𝜅
+
𝜅
− 𝑤𝑄 𝑥, 𝑦 𝑑𝑥𝑑𝑦
𝜕𝑥
𝜕𝑥
𝜕𝑦
𝜕𝑦
Ω
−
𝑤
Γ
𝜕𝑇
𝜕𝑇
𝜅
𝑛 + 𝜅
𝑛 𝑑𝑠 = 0
𝜕𝑥 𝑥
𝜕𝑦 𝑦
EBC: Specify T(x,y) on G
T
T
NBC: Specify x nx y n y on G
T T
where qn ( s) i j nxi n y j is the normal
x y
outward flux on the boundary G at the segment ds.
FEM Implementation of 2-D Heat
Conduction – Shape Functions
Step 1: Discretization – linear triangular element (T3)
T1
T T11 T22 T33
Derivation of linear triangular shape functions:
T2
1 c0 c1 x c2 y
Let
T3
Interpolation properties
i 1 at ith node
i 0 at other nodes
1 1 x
Same
2
1
1 x1
y 1 x2
1 x3
c0 c1 x1 c2 y1 1
c0 c1 x2 c2 y2 0
c0 c1 x3 c2 y3 0
y1
y2
y3
x3 y1 x1 y3
y
y3 y1
2 Ae
x x
1 3
x
c0 1 x1
c 1 x
2
1
c 1 x
3
2
1
1
x2 y3 x3 y2
1
x
y
0
y
y
2
3
2
A
e
0
x3 x2
3
1
x1 y2 x2 y1
y
y1 y2
2 Ae
x x
2 1
x
y1
y2
y3
1
1
0
0
FEM Implementation of 2-D Heat
Conduction – Shape Functions
linear triangular element – local (area) coordinates
T1
𝜙1 =
A2
T2
1
𝑥 𝑦
2𝐴𝑒
𝑥2 𝑦3 − 𝑥3 𝑦2
𝐴1
𝑦2 − 𝑦3
=
=𝜉
𝐴
𝑒
𝑥3 − 𝑥2
1
𝑥 𝑦
2𝐴𝑒
𝑥3 𝑦1 − 𝑥1 𝑦3
𝐴2
𝑦3 − 𝑦1
=
=𝜂
𝐴
𝑒
𝑥1 − 𝑥3
1
𝑥 𝑦
2𝐴𝑒
𝑥1 𝑦2 − 𝑥2 𝑦1
𝐴3
𝑦1 − 𝑦2
=
=1−𝜉−𝜂 =𝜁
𝐴
𝑒
𝑥2 − 𝑥1
A3
A1
𝜙2 =
T3
T1
𝜉=1
𝜂=1
𝜙3 =
T2
T3
𝜉=0
𝜂=0
1
3
2
FEM Implementation of 2-D Heat
Conduction – Shape Functions
quadratic triangular element (T6) – local (area) coordinates
T1
T4
T6
𝜙1 = 𝜉 2𝜉 − 1
T2
𝜙2 = 𝜂 2𝜂 − 1
T3
T5
𝜙3 = 𝜁 2𝜁 − 1
𝜙4 = 4𝜉𝜂
T1
𝜂=1
𝜉=1
T4
T6
𝜉 = 0.5
T2
T5
T3
𝜙5 = 4𝜂𝜁
𝜙6 = 4𝜉𝜁
𝜉=0
𝜂=0
Serendipity Family – nodes are placed on the boundary
for triangular elements, incomplete beyond quadratic
Interpolation Function - Requirements
• Interpolation condition
• Take a unit value at node i, and is zero at all other nodes
• Local support condition
• i is zero at an edge that doesn’t contain node i.
• Interelement compatibility condition
• Satisfies continuity condition between adjacent elements
over any element boundary that includes node i
• Completeness condition
• The interpolation is able to represent exactly any
displacement field which is polynomial in x and y with the
order of the interpolation function
Formulation of 2-D 4-Node Rectangular Element –
Bi-linear Element (Q4)
Let u( , ) 1u1 2u2 3u3 4u4
4
3
1
1
4
1
𝜙2 =
4
1
𝜙3 =
4
1
𝜙4 =
4
𝜙1 =
2
1−𝜉 1−𝜂
1+𝜉 1−𝜂
1+𝜉 1+𝜂
1−𝜉 1+𝜂
Note: The local node numbers should be arranged in a counter-clockwise sense. Otherwise, the area
Of the element would be negative and the stiffness matrix can not be formed.
1
2
3
4
Formulation of 2-D 4-Node Rectangular Element –
Bi-linear Element (Q4)
Physical domain (physical element)
4
3
Reference domain (master element)
4
3
y
1
x
2
1
𝑥 = 𝜙1 𝑥1 + 𝜙2 𝑥2 + 𝜙3 𝑥3 + 𝜙4 𝑥4
𝑦 = 𝜙1 𝑦1 + 𝜙2 𝑦2 + 𝜙3 𝑦3 + 𝜙4 𝑦4
2
FEM Implementation of 2-D Heat
Conduction – Element Equation
• Weak Form of 2-D Model Problem ----𝜕𝑤
𝜕𝑇
𝜕𝑤
𝜕𝑇
𝜅
+
𝜅
− 𝑤𝑄 𝑥, 𝑦 𝑑𝑥𝑑𝑦 +
𝜕𝑥
𝜕𝑥
𝜕𝑦
𝜕𝑦
Ω𝑒
𝑛
Assume approximation:
𝑇=
𝑤𝑞𝑛 𝑑𝑠 = 0
Γ𝑒
𝑇𝑗 𝜙𝑗
𝑗=1
and let w(x,y)=i(x,y) as before, then
Ω𝑒
𝜕𝜙𝑖 𝜕
𝜅
𝜕𝑥 𝜕𝑥
𝑛
𝑇𝑗 𝜙𝑗
𝑗=1
𝜕𝜙𝑖 𝜕
+
𝜅
𝜕𝑦 𝜕𝑦
𝑛
𝑇𝑗 𝜙𝑗 − 𝜙𝑖 𝑄 𝑥, 𝑦
𝑑𝑥𝑑𝑦
𝑗=1
𝑛
+
𝜙𝑖 𝑞𝑛 𝑑𝑠 = 0
𝐾𝑖𝑗 𝑇𝑗 =
𝑗=1
Γ𝑒
where
i j
i j
K ij
dxdy
x x
y y
We
𝜙𝑖 𝑄 𝑥, 𝑦 𝑑𝑥𝑑𝑦 −
Ω𝑒
𝜙𝑖 𝑞𝑛 𝑑𝑠
Γ𝑒
Linear Triangular Element Equation
𝑛
𝐾𝑖𝑗 𝑇𝑗 =
𝑗=1
1
Ω𝑒
x2 y3 x3 y2
y
A1
y2 y3
2 Ae
x x Ae
3 2
1 x
2
3
𝜙𝑖 𝑄 𝑥, 𝑦 𝑑𝑥𝑑𝑦 −
x3 y1 x1 y3
y
A2
y3 y1
2 Ae
x x Ae
1 3
1 x
1
x1 y2 x2 y1
y
A3
y1 y2
2 Ae
x x Ae
2 1
x
𝜙𝑖 𝑞𝑛 𝑑𝑠
Γ𝑒
l23 l23
K l23 l31
4 Ae
l23 l12
l23 l31 l23 l12
l31 l31 l31 l12
l31 l12 l12 l12
Q1 q1
F Q2 q2
Q q
3 3
where 𝒍𝒊𝒋 is the length vector from the ith node
to the jth node.
Assembly of Stiffness Matrices
𝑛𝑒
𝐹𝑖
𝑒
=
𝜙𝑖 𝑄 𝑥, 𝑦 𝑑𝑥𝑑𝑦 −
Ω𝑒
𝑈1 = 𝑇1
1
, 𝑈2 = 𝑇2
𝜙𝑖 𝑞𝑛 𝑑𝑠 =
+ 𝑇1
2
, 𝑈3 = 𝑇3
𝑒
𝑇𝑗
𝑒
𝑗=1
Γ𝑒
1
𝐾𝑖𝑗
1
+ 𝑇4
2
, 𝑈4 = 𝑇2
2
, 𝑈5 = 𝑇3
2
,
Imposing Boundary Conditions
The meaning of qi:
q1(1)
3
3
1
1
(1) (1)
q
n 1 ds
G1
1
(1)
h12
2
(1) (1)
q
n 2 ds
G1
(1)
h12
3
1
1
2
(1)
h23
(1)
h12
qn(1)2(1) ds qn(1)2(1) ds qn(1)2(1) ds
(1)
h23
(1)
h31
(1)
h23
(1) (1)
q
n 3 ds
G1
2
(1)
h31
qn(1)2(1) ds qn(1)2(1) ds
q3(1)
1
(1)
h23
(1)
h31
q2(1)
3
(1)
h12
qn(1)1(1) ds qn(1)1(1) ds
2
qn(1)1(1) ds qn(1)1(1) ds qn(1)1(1) ds
(1)
h12
qn(1)3(1) ds qn(1)3(1) ds qn(1)3(1) ds
qn(1)3(1) ds qn(1)3(1) ds
(1)
h31
(1)
h23
(1)
h31
Imposing Boundary Conditions
Consider
(1)
2
q
q
ds q ds
q1(2)
q
ds q ds
q4(2)
(1) (1)
n
2
(1)
h12
(1)
3
q2 q q
q3 q3(1) q4(2)
(1)
2
(2)
1
(1) (1)
n
2
q
(1)
h23
(1) (1)
n
3
Equilibrium of flux:
(1)
h23
qn(2)
qn(2)1(2) ds
qn(2)4(2) ds
(2)
h41
(2)
h34
(1)
h31
qn(1)
qn(2)4(2) ds
(2)
h12
(1)
h23
(1) (1)
n
3
qn(2)1(2) ds
(2)
h41
( 2)
h41
FEM implementation:
qn(1)2(1) ds
(1)
h23
qn(2)1(2) ds;
( 2)
h41
q2
(1)
h12
qn(1)2(1) ds
qn(1)3(1) ds
(1)
h23
( 2)
h12
qn(2)1(2) ds
qn(2)4(2) ds
( 2)
h41
q3
(1)
h31
qn(1)3(1) ds
( 2)
h34
qn(2)4(2) ds
Calculating the q Vector
Example:
qn 0
T 293K
qn 1
2-D Steady-State Heat Conduction - Example
A
D
AB:
qn 0
CD: convection
h 50W
DA and BC: T 180 C
o
0.6 m
C
B
0.4 m
y
x
m C
2 o
T 25o C