Lecture 6

Report
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  T11  T22  T33
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

similar documents