Spring Element The first step in FEM is the derivation of the element stiffness matrix ke. This is illustrated below by using a simple spring element x F1
2 F2
1
u1
k
u2
Two nodes: 1,2 Nodal displacement: u1, u2 (m, mm) Nodal forces: F1, F2 (N) Spring constant (Stiffness) k (N/m, N/mm) 7
Spring Element We have:
F =kΔ
with Δ =u2 - u1
(1)
Consider the equilibrium of forces for the spring. At node 1, we have:
and at node 2, In matrix form
or
F1 = − F = −k(u2 − u1) = ku1 − ku2
(2)
F2 = F = k(u2 − u1) = −ku1 + ku2
(3)
⎡ k − k⎤ ⎧u1 ⎫ ⎧ F1 ⎫ ⎢− k k ⎥ ⎨u ⎬ = ⎨F ⎬ ⎣ ⎦ ⎩ 2⎭ ⎩ 2⎭
(4)
ku =F
(5)
Where k =(element) stiffness matrix u =(element nodal) displacement vector F = (element nodal) force vector
8
Bar Element Consider a uniform prismatic bar:
L A E
length cross-sectional area elastic modulus 9
Stiffness matrix – Direct method Assume that the displacement in the bar is δ, then we have:
and
δ = FL / AE
(6)
F = ( AE / L)δ
(7)
For the respective cases in the above figure, we have: F11 = F21 =
AE u1 L
and
F12 = F22 =
AE u2 L
(8)
where Fij is the force at node i (i = 1,2) associated with displacement of node j (j = 1,2)
10
Stiffness matrix – Direct method Written the above expressions in matrix format:
⎡ F11 − F12 ⎤ ⎧1⎫ ⎧ F1 ⎫ ⎨ ⎬=⎨ ⎬ ⎢− F ⎥ ⎣ 21 F22 ⎦ ⎩1⎭ ⎩F2 ⎭
or
AE ⎡ 1 − 1⎤ ⎧u1 ⎫ ⎧ F1 ⎫ ⎨ ⎬=⎨ ⎬ ⎢ ⎥ L ⎣− 1 1 ⎦ ⎩u2 ⎭ ⎩ F2 ⎭
(9)
where F1 and F2 are the resultant forces applied to the bar at nodes 1 and 2 From the above matrix, we can see that the bar is acting like a spring in this case, the element stiffness matrix for the bar is:
⎡ k − k⎤ AE ⎡ 1 − 1⎤ = k= ⎢ ⎥ ⎢− 1 1 ⎥ − k k L ⎣ ⎦ ⎣ ⎦
(10) 11
Stiffness matrix – A Formal Procedure • We derive the same stiffness matrix for the bar using a formal procedure which can be applied to more complicated situations • For different types of element, the variation of the displacement between the nodes must be assumed first (what kind of behaviour can be expected). This is called thedisplacement function and is usually a polynomial whose order depends on the number of degree of freedom (d.o.f.) in the element. Degree of freedom is number of components of the displacement vector at a node. 12
Stiffness matrix – A Formal Procedure Bar elements are used to model truss structures and any other structures where axial effects predominate. Bar element consists of 2 d.o.f. per element and hence the displacement function can be assumed as: u=
1
+ β 2x
(11)
The FEM treats the nodal displacements as variables of an interpolation function, usually a polynomial, to give an analytical expression for displacement at any point inside the element. 13
Stiffness matrix – A Formal Procedure Eq. (11) can be written by expressing the βi in terms of nodal displacements u1 and u2: at x =0 at x =L
u =u1 =β1 u =u2 =β1+β2L
replace β1 and β2 into eq. (11), we have:
x x u = (1− )u1 + u2 L L
(12)
Define two linear shape functions as follows: N1(ξ ) = 1− ξ
N2(ξ ) = ξ
(13)
where
x ξ = L
0≤ ξ ≤ 1
14
Stiffness matrix – A Formal Procedure Each shape function Ni describes how u varies with x when the corresponding d.o.f. ui is unity while the other is zero. Re-write eq. (12) in matrix form: L − x x⎤ ⎧u1 ⎫ ⎡ u= ⎢ or ⎨ ⎬ u =Nd ⎥ u L L ⎣ ⎦⎩ 2 ⎭ where N is theshape function matrix d is the vector of element nodal d.o.f.
(14)
15
Stiffness matrix – A Formal Procedure Axial strain εx is given by ε =
du ⎡ d ⎤ = ⎢ N ⎥ d = Bd dx ⎣ dx ⎦
(15)
whereB is the element strain-displacement matrix ⎡ 1 1⎤ B = ⎢− ⎣ L L ⎥⎦
(16)
Stress can be written as
σ =Eε =EBd
(17) 16
Stiffness matrix – A Formal Procedure Consider thestrain energy stored in the bar U=
1 T 1 T T σ ε dV = (d B EBd)dV 2V 2V
∫
∫
⎤ 1 T ⎡ T = d ⎢ (B EB)dV ⎥ d 2 ⎣V ⎦
(18)
∫
Thework done by the two nodal forces is W=
1 1 1 F1u1 + F2u2 = dT F 2 2 2
(19)
For conservative system, we have U =W
(20)
17
Stiffness matrix – A Formal Procedure which gives ⎤ 1 T ⎡ T 1 d ⎢ (B EB)dV ⎥ d = dT F 2 ⎣V 2 ⎦
∫
(21)
Then we have ⎡ T ⎤ ⎢ ∫ (B EB)dV ⎥ d = F ⎣V ⎦
where
∫
(22)
k = (B T EB)dV
or
kd =F
(23)
(24)
V
is theelement stiffness matrix Equation (24) is a general result which can be used for the construction of other types of elements.
18
Stiffness matrix – A Formal Procedure Now we use equation (24) to re-calculate the element stiffness matrix for the bar element:
⎧− 1/ L⎫ AE ⎡ 1 − 1⎤ k= ∫ ⎨ ⎬E[− 1/ L 1/ L]Adx = ⎢− 1 1 ⎥ 0 1 / L L ⎩ ⎭ ⎣ ⎦ L
(25)
which is the same as we derived using the direct method Thestrain energy can be written as 1 T U = d kd 2
(26) 19
Example: Find the stresses in the two bar assembly which is loaded with force P, and constrained at the two ends, as shown in the figure.
20
Notes to example • In this case, the calculated stresses in elements 1 and 2 are exact within the linear theory for for 1D bar structures. It will not help if we further divide element 1 or 2 into smaller finite elements. • For tapered bars, averaged values of the crosssectional areas should be used for the elements. • We need to find the displacements first in order to find the stresses, since we are using the displacement based FEM. 21
Beam Element The beam element is used to model beams or frames where flexural effects (shear forces and bending moments) dominate. Beam element consists of 4 DOF per element and a cubic variation in displacement has to be assumed in the form: v = β 1 + β 2 x + β 3x2 + β 4 x3 (27) I =moment of inertia of the crosssectional area E = elastic modulus v =v(x) lateral displacement dv rotation about the z-axis θ = dx F = shear force M =bending moment about z-axis
22
Beam Element Axial effect can be added on if necessary Similar to bar element, eq. (27) can be written by expressing the βi in terms of nodal d.o.f. with θz=dv/dx=β2+2β3x+3β4x2 for example: x =0, β1 =v1, β2 =θz1 We can then derive the shape functions for beam element. N1( x) = 1− 3x2 / L2 + 2x3 / L3 N2 ( x) = x − 2x2 / L + x3 / L2 N3( x) = 3x / L − 2x / L 2
2
3
N4 ( x) = x2 / L + x3 / L2
3
(28) 23
Beam Element Then, we can represent the displacement as v = [N1
N2
N3
⎧ v1 ⎫ ⎪θ ⎪ ⎪ z1 ⎪ N4 ]⎨ ⎬ = Nd ⎪ v2 ⎪ ⎪⎩θ z2 ⎪⎭
(29)
Curvature of the beam element is d 2v ⎡ d 2 ⎤ = ⎢ 2 N ⎥ d = Bd 2 dx ⎣ dx ⎦
(30)
where strain-displacement matrix B is given by 4 6x 6 12x 2 6x⎤ ⎡ 6 12x + − + − − + 2⎥ 2 3 2 2 3 L L L L ⎦ L L L ⎣ L
B = ⎢−
(31)
24
Beam Element Strain energy stored in the beam element is
⎤ 1 T 1 T ⎡ T U= σ ε dV = d ⎢ B EIBdx⎥ d 2V 2 ⎣0 ⎦ L
∫
∫
(32)
Can you prove this?? Based on eq (32) the stiffness matrix for the beam element is: L
∫
k = B T EIBdx
(33)
0
25
Beam Element Applying the result in (31) and carrying out the integration, the expression of k is as follows ⎡ 12EI / L3 6EI / L2 − 12EI / L3 6EI / L2 ⎤ ⎢ ⎥ 2 2 6 EI / L 4 EI / L 6 EI / L 2 EI / L − ⎥ k= ⎢ ⎢− 12EI / L3 − 6EI / L2 12EI / L3 − 6EI / L2 ⎥ ⎢ ⎥ 2 2 6 EI / L 2 EI / L 6 EI / L 4 EI / L − ⎣⎢ ⎦⎥
Stress
(34)
σx =My/I
Bending moment M is computed from curvature d2v/dx2, which in turn depends on nodal d.o.f. d. d 2v M = EI 2 = EIBd dx
(35)
26
2D Beam Element Combining the axial stiffness − AE / L 0 0 0 0 ⎡ AE / L ⎤ u1 ⎢ 0 3 2 3 2 ⎥ EI / L EI / L EI / L EI / L − 12 6 0 12 6 ⎢ ⎥ v1 ⎢ 0 − 6EI / L2 2EI / L ⎥ θ z1 6EI / L2 4EI / L 0 k= ⎢ ⎥ AE / L AE / L − 0 0 0 0 ⎢ ⎥ u2 ⎢ 0 − 12EI / L3 − 6EI / L2 0 12EI / L3 − 6EI / L2 ⎥ v2 ⎢ ⎥ 2 2 EI / L EI / L EI / L EI / L − 0 6 2 0 6 4 ⎣⎢ ⎦⎥θ z2
27
Boundary conditions • In order to prevent the finite element model from moving freely through space, each of the possible degrees of freedom must be constrained somewhere on the model. This can be done by imposing boundary conditions at appropriate nodal points on the model. • In general, displacement boundary conditions simulate the actual supports of the structure. Same as general structural analysis, there are three basic types of supports: simply supported, fixed and roller supports. 28
Boundary Conditions • The application of boundary conditions depends on the type of structure being analysed and also on the finite element program being used. • If the program supports 2D elements, then the application of boundary conditions in the third dimension for a 2D analysis is unnecessary. On the other hand, programs that use 3D elements for 2D problems require that all displacements in the third dimension are zeroed. 29
