05-Bar Element Duzeltilmis

ME 520
Fundamentals of Finite Element

Analysis
5-Bar Element
Dr. Ahmet Zafer enalp

e-mail: [email protected]
Mechanical Engineering Department

Gebze Technical University
5-Bar Element
Bar Element
Bar (truss) structures:
ME 520
Mechanical Engineering Department, GTU
5-Bar Element
Bar Element
Consider a uniform prismatic bar:
FE Model
L: length
A: cross-sectional area
E: elastic modulus
u=u(x): displacement
=(x): strain
=(x): stress
Strain-displacement relation:
Stress-strain relation:
ME 520
Stiffness Matrix - Direct Method
5-Bar Element
Assuming that the displacement u is varying linearly along the axis of

the bar, i.e.,
: elongation, F: force in bar

k: stiffness of the bar
The bar is acting like a spring in this case and we conclude that
element stiffness matrix is:
ME 520
Stiffness Matrix - A Formal Approach
5-Bar Element
We derive the same stiffness matrix for the bar using a formal
approach which can be applied to many other more complicated
situations:
Define two linear shape functions as follows:
B: Element strain-displacement matrix
ME 520
5-Bar Element
Stress can be written as:

Consider the strain energy stored in the bar:
The work done by the two nodal forces is:
For conservative system, we state that:

which gives:
ME 520
5-Bar Element
k is the element stiffness matrix. The above expression is a general

result which can be used for the construction of other types of
elements. This expression can also be derived using other more
rigorous approaches, such as the Principle of Minimum Potential
Energy, or the Galerkins Method.
Now, we evaluate for the bar element
which is the same as we derived using the direct method.

Strain energy in the element can be written as:
ME 520
Degree of Freedom (dof)
5-Bar Element
Number of components of the displacement vector at a node.

For 1-D bar element: one dof at each node.
Physical Meaning of the Coefficients in k:
The jth column of k (here j = 1 or 2) represents the forces
applied to the bar to maintain a deformed shape with unit
displacement at node j and zero displacement at the other node.
ME 520
5-Bar Element
Example 1
FE Model
1
Find;
(a) the global stiffness matrix
(b) displacements of nodes
(c) the reaction forces
(d) stresses at the elements
Solution :
Connectivity table:
ME 520
E
#
N1
N2
5-Bar Element
Example 1
Boundary conditions :
Displacement boundary conditions
u1 0,:u 2 0, u 3 0
F1 :0, F2 P, F3 0
Force boundary conditions
a) Element Stiffness Matrices :
ME 520
10
5-Bar Element
Example 1
Construction of global stiffness matrix :
Equilibrium (FE)
equation for
the whole system is;
ME 520
11
5-Bar Element
Example 1
u1 0, u 2 0, u 3 0 F1 0, F2 P, F3 0
b) Applying boundary conditions;
ME 520
12
5-Bar Element
Example 1
For the data given;

P 5x10 4 N, E 2x10 4 N / mm 2
A 200 mm 2 , L 125 mm
0
u1
u 2 0.5208 mm
u
0
3
ME 520
13
5-Bar Element
Example 1
c) From the 1st and 3rd row of the finite element equation, the reaction
forces can be calculated as;
F1 -33333.3 N
F2 -16666.7 N
d) Stresses in the elements:
1 83.33 MPa
2 83.33 MPa
ME 520
14
5-Bar Element
Notes
In this case, the calculated stresses in elements 1 and 2 are exact

within the linear theory for 1-D 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 cross-sectional 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.
ME 520
15
Solution procedure with matlab
5-Bar Element
As linear bar element has 2 degrees of freedom (1 at each

node) for a structure with n nodes, the global stiffness matrix K
will be of size nxn.
The global stiffness matrix K is obtained by making calls to the
Matlab function LinearBarAssemble which is written for this
purpose.
Once the global stiffness matrix; K is obtained we have the
following
equation;
K U structure
F
At this step boundary conditions are applied manually to the

vectors U and F.
Then the matrix equation is solved by partioning and Gaussion
elimination.
ME 520
Finally
once the unknown

andGTU
and reactions16are
Mechanical displacements
Engineering Department,
found, the force is obtained for each element as follows:
Solution procedure with matlab
5-Bar Element
f k u
where f is the 2x1 element force vector and u is the 2x1
element displacement vector.
The element stress is obtained by using local stiffness matrix
and local displacement vector.
ME 520
17
Matlab functions used
5-Bar Element
LinearBarElementStiffness(E,A,L)
This function returns the element stiffness matrix for a linear bar
with modulus of elasticity E, cross-sectional area A, and length L.
The size of the element stiffness matrix is 2 x 2.
Function contents:
function y = LinearBarElementStiffness(E,A,L)
%LinearBarElementStiffness
This function returns the element
%
stiffness matrix for a linear bar with
%
modulus of elasticity E, cross-sectional
%
area A, and length L. The size of the
%
element stiffness matrix is 2 x 2.
y = [E*A/L -E*A/L ; -E*A/L E*A/L];
ME 520
18
5-Bar Element
LinearBarAssemble(K,k,i,j)
This function assembles the element stiffness matrix k of the linear
bar with nodes i and j into the global stiffness matrix K. This
function returns the global stiffness
matrix K after the
element stiffness matrix k is assembled.
Function contents:
function y = LinearBarAssemble(K,k,i,j)
%LinearBarAssemble
This function assembles the element stiffness
%
matrix k of the linear bar with nodes i and j
%
into the global stiffness matrix K.
%
This function returns the global stiffness
%
matrix K after the element stiffness matrix
%
k is assembled.
K(i,i) = K(i,i) + k(1,1);
K(i,j) = K(i,j) + k(1,2);
K(j,i) = K(j,i) + k(2,1);
K(j,j) = K(j,j) + k(2,2);
y = K;
ME 520
19
5-Bar Element
LinearBarElementForces(k,u)
This function returns the element nodal force vector given the
element stiffness
matrix k and the element nodal displacement vector u.
Function contents:
function y = LinearBarElementForces(k,u)
%LinearBarElementForces
This function returns the element nodal
%
force vector given the element stiffness
%
matrix k and the element nodal displacement
%
vector u.
y = k * u;
ME 520
20
5-Bar Element
LinearBarElementStressess_az(k, u, L, E)
This function returns the element nodal stress vector given the
element stiffness
matrix k, the element nodal
displacement vector u, the length of the element and E
Function contents:
function y = LinearBarElementStresses_az(k, u, L, E)
%LinearBarElementStresses
This function returns the element nodal
%
stress vector given the element stiffness
%
matrix k, the element nodal displacement
%
vector u, the length of the element and E
L_m=[-1/L 1/L];
y = E*L_m*u;
ME 520
21
5-Bar Element
Solution of Example 1
with Matlab
Solution:
Use the 7 steps to solve the problem using linear bar element.
Step 1-Discretizing the domain:
This problem is already discretized. The domain is subdivided into
2 elements and 3 nodes.
The units used in Matlab calculations are N and milimeter. The
element connectivity is:
ME 520
E
#
N1
N2
22
with Matlab
5-Bar Element
P 5x10 4 N, E 2x10 4 N / mm 2
A 200 mm 2 , L 125 mm
Step 2-Copying relevant files and starting Matlab

Create a directory
Copy
LinearBarElementStiffness.m
LinearBarAssemble.m
LinearBarElementForces.m
LinearBarElementStresses_az.m
files under the created directory
Open Matlab;
Open Set Path command and by using Add Folder command add the
current directory.
ME 520
23
with Matlab
5-Bar Element
P 5x10 4 N, E 2x10 4 N / mm 2
A 200 mm 2 , L 125 mm
Start solving the problem in Command Window:

>>clearvars
>>clc
Enter the data
>>E=2e4
>>A(1)=400
>>A(2)=200
>>L=125
Step 3-Writing the element stiffness matrices:
>>k1=LinearBarElementStiffness(E,A(1),L)
k1 =
64000
-64000
-64000
64000
ME 520
24
with Matlab
5-Bar Element
P 5x10 4 N, E 2x10 4 N / mm 2
A 200 mm 2 , L 125 mm
k2 =
32000
-32000
-32000
32000
Step 4-Assembling the global stiffness matrix:

Since the structure has 3 nodes, the size of the global stiffness matrix is
3x3.
Therefore to obtain K we first set up a zero matrix of size 3x3 then make
2 calls to the Matlab function LinearBarAssemble since we have 2 bar
elements in the system.
Each call to the function will assemble one element. The following are
the Matlab commands:
ME 520
25
5-Bar Element
with Matlab
P 5x10 4 N, E 2x10 4 N / mm 2
A 200 mm 2 , L 125 mm
>>K=zeros(3,3)
K=
0
0
0
0
0
0
0
0
0
>>K=LinearBarAssemble(K,k1,1,2)
K=
64000
-64000
0
-64000
64000
0
0
0
0
K=
64000
-64000
0
-64000
96000
-32000
0
-32000
32000
ME 520
26
5-Bar Element
with Matlab
P 5x10 4 N, E 2x10 4 N / mm 2
A 200 mm 2 , L 125 mm
Step 5-Applying the boundary conditions:

Finite element equation for the problem is;
0
64000 - 64000
- 64000 96000 - 32000
0
- 32000 32000
u1
F1

u 2 F2
u
F
3
3
The boundary conditions for the problem are;
u1 0, u 2 0, u 3 0
F1 0, F2 P, F3 0
Inserting the above conditions into finite element equation

0
64000 - 64000
- 64000 96000 - 32000
0
- 32000 32000
ME 520
0
F1
u 2 50000
0
F
3
27
5-Bar Element
with Matlab
Step 6-Solving the equations:

Solving the above system of equations will be performed by
partitioning (manually) and Gaussian elimination (with Matlab) First
we partition the above equation;
0
64000 - 64000
- 64000 96000 - 32000
0
- 32000 32000
0
F1
u 2 50000
0
F
3
96000u 2 50000
>>u2=50000/96
000
u 0.5208
2
Step 7-Post-processing:
In this step we obtain the reactions at nodes 1 and 3 and the
stresses at each element using Matlab as follows.
First we set up the global nodal displacement vector U, then we
calculate the nodal force vector F.
ME 520
28
with Matlab
5-Bar Element
>>U=[0; u2; 0]
U=
0
0.5208
0
>>F=K*U
F=
1.0e+04 *
-3.3333
5.0000
-1.6667
ME 520
29
with Matlab
5-Bar Element
P 5x10 4 N, E 2x10 4 N / mm 2
A 200 mm 2 , L 125 mm
>>u_el_1=[U(1); U(2)]
>>u_el_2=[U(2); U(3)]
>>stress1=LinearBarElementStresses_az(k1, u_el_1, L, E)
>>stress2=LinearBarElementStresses_az(k2, u_el_2, L, E)
stress1 =
83.3333
:element stresses
stress2 =
-83.3333
ME 520
30
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
node 3 is displayed to the right
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
Create a directory
Copy
LinearBarAssemble.m
Open Matlab;
Open Set Path command and by using Add Folder command add the
current directory.
ME 520
31
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m

>>clearvars
>>clc
Enter the data
>>E=210e6
>>A=0.003
>>L1=1.5
>>L2=1
Calculate stiffness matrices
>>k1=LinearBarElementStiffness(E,A,L1)
>>k2=LinearBarElementStiffness(E,A,L2)
ME 520
32
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
Assemble the global stiffness matrix;

>>K=zeros(3,3)
K is calculated as:
K=
420000
-420000
0
-420000
1050000
-630000
0
-630000
630000
ME 520
33
5-Bar Element
with Matlab
P 10 kN, E 210 GPa
u1 0, u 2 0, u 3 0.002
A 0.003 m 2
F1 0, F2 P, F3 0
by 0.002 m
FE Equation is:
0
420000 - 420000
- 420000 1050000 - 630000
0
- 630000 630000
u1
F1

u2 F2
u
F
3
3
Apply BCs:u1 0, u 2 0, u 3 0.002

0
420000 - 420000
- 420000 1050000 - 630000
0
- 630000 630000
ME 520
F1 0, F2 P, F3 0
0
F1
u2 10
0.002
F
34
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
Solving the above system of equations will be performed by partitioning

(manually) and Gaussian elimination (with Matlab).
First we partition the above equation by extracting the submatrix in row
2 and column 2 which turns to be a 1x1 matrix.
Because of the applied displacement of 0.002 m at node 3, we need to
extract the submatrix in row 2 and column 3 which also turns out to be
a 1x1 matrix.
Therefore we obtain;
1050000 u 2 630000 0.002 10
ME 520
35
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
The solution of the above system is obtained using Matlab as follows.

Note that the backslash operator \ is used for Gaussian elimination.
>>k=K(2,2)
>>k0=K(2,3)
>>u0=0.002
>>f=[-10]
>>f0=f-k0*u0
>>u=k\f0
result is;
u=
0.0012
ME 520
36
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
In this step, we obtain the reactions at nodes 1 and 3, and the stress in
each bar using Matlab as follows.
First we set up the global displacement vector; U, then we calculate the
golbal force vector F.
>>U=[0; u; u0]
>>F=K*U
results;
F=
-500.0000
-10.0000
510.0000
ME 520
37
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
Next step is to set up element nodal displacement vectors u 1 and u2.

Then we calculate the element force vectors f 1 and f2 by making calls to
the Matlab function LinearBarElementForces.
>>u1=[0; U(2)]
>>f1=LinearBarElementForces(k1,u1)
result is;
f1 =
-500.0000
500.0000
ME 520
38
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
>>u2=[U(2); U(3)]
>>f2=LinearBarElementForces(k2,u2)
result is;
f2 =
-510.0000
510.0000
Finally, we call LinearBarElementStresses_az to calculate stresses.
>>sigma1=LinearBarElementStresses_az(k1, u1, L1, E)
results;
sigma1 =
1.6667e+05
ME 520
39
with Matlab
P 10 kN, E 210 GPa
A 0.003 m 2
5-Bar Element
u1 0, u 2 0, u 3 0.002
F1 0, F2 P, F3 0
by 0.002 m
>>sigma2=LinearBarElementStresses_az(k2, u2, L2, E)

results;
sigma2 =
1.7000e+05
ME 520
40
5-Bar Element
Example 3
Given;
Find;
The support reaction forces at the two ends of the bar
Solution :
Connectivity table:
ME 520
E#
N1
N2
41
5-Bar Element
Example 3
We first check to see if or not the contact of the bar with the wall on
the right will occur. To do this, we imagine the wall on the right is
removed and calculate the displacement at the right end,
contact occurs
BCs; u1 0, u 2 0, u 3
F1 0, F2 P, F3 0
The global FE equation is found to be,
ME 520
42
5-Bar Element
Example 3
u1 0, u 2 0, u 3
Applying boundary conditions;
ME 520
F1 0, F2 P, F3 0
43
5-Bar Element
Example: 3
Reaction forces:
ME 520
44
5-Bar Element
with Matlab
E
#
N1
N
2
contact occurs
BCs;
u1 0, u 2 0, u 3
F1 0, F2 P, F3 0
ME 520
45
with Matlab
5-Bar Element
E
#
N1
N
2
u1 0, u 2 0, u 3
F1 0, F2 P, F3 0
Create a directory
Copy
LinearBarAssemble.m
Open Matlab;
Open Set Path command and by using Add Folder
command add the current directory.
ME 520
46
with Matlab
5-Bar Element
E
#
N1
N
2
u1 0, u 2 0, u 3
F1 0, F2 P, F3 0
>>clearvars
>>clc
Enter the data
>>P=6e4
>>E=2e4
>>A=250
>>L=150
>>DELTA=1.2
Calculate stiffness matrices
>>k1=LinearBarElementStiffness(E,A,L)
>>k2=LinearBarElementStiffness(E,A,L)
ME 520
47
with Matlab
5-Bar Element
E
#
N1
N
2
u1 0, u 2 0, u 3
Assemble the global stiffness matrix;
>>K=zeros(3,3)
K is calculated as:
F1 0, F2 P, F3 0
K=
1.0e+04 *
3.3333 -3.3333
0
-3.3333 6.6667 -3.3333
0 -3.3333 3.3333
ME 520
48
5-Bar Element
with Matlab
E
#
N1
N
2
u1 0, u 2 0, u 3
F1 0, F2 P, F3 0
FE Equation is:
0
33333 - 33333
- 33333 66667 - 33333
0
- 33333 33333
u1
F1

u 2 F2
u
F
3
3
Apply BCs:u1 0, u 2 0, u 3
0
33333 - 33333
- 33333 66667 - 33333
0
- 33333 33333
ME 520
F1 0, F2 P, F3 0
0
F1
u 2 60000
1 .2
F
3
49
with Matlab
5-Bar Element
E
#
N1
N
2
u1 0, u 2 0, u 3
F1 0, F2 P, F3 0
66667 u 2 333331.2 60000

>>k=K(2,2)
>>k0=K(2,3)
>>u0=1.2
>>f=[60000]
>>f0=f-k0*u0
>>u=k\f0
result is;
u=
1.5000
ME 520
50
with Matlab
5-Bar Element
E
#
N1
N
2
u1 0, u 2 0, u 3
F1 0, F2 P, F3 0
>>U=[0; u; u0]
>>F=K*U
results;
F=
-50000
60000
-10000
ME 520
51
with Matlab
5-Bar Element
E
#
N1
N
2
u1 0, u 2 0, u 3
>>u1=[0; U(2)]
>>f1=LinearBarElementForces(k1,u1
)
results;
f1 =
F1 0, F2 P, F3 0
-50000
50000
>>u2=[U(2); U(3)]
>>f2=LinearBarElementForces(k2,u2
)
results;
f2 =
ME 520
10000
52
with Matlab
5-Bar Element
E
#
N1
N
2
u1 0, u 2 0, u 3
>>sigma1=LinearBarElementStresses_az(
k1, u1, L, E)
results;
sigma1 =
F1 0, F2 P, F3 0
200
>>sigma2=LinearBarElementStresses_az(
k2, u2, L, E)
results;
sigma2 =
-40
ME 520
53
5-Bar Element
Distributed Load
Uniformly distributed axial load q (N/mm, N/m, lb/in) can be converted

to two equivalent nodal forces of magnitude qL/2.
We verify this by considering the work done by the load q,
ME 520
54
5-Bar Element
Distributed Load
Thus, from the U=W concept for the element, we have:
The new nodal force vector is:
ME 520
55
5-Bar Element
Distributed Load
ME 520
56

05-Bar Element Duzeltilmis

Uploaded by

Copyright:

Available Formats

05-Bar Element Duzeltilmis

Uploaded by

Document Information

Original Description:

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

05-Bar Element Duzeltilmis

Uploaded by

Copyright:

Available Formats

ME 520

Fundamentals of Finite Element

Dr. Ahmet Zafer enalp

Mechanical Engineering Department

Bar (truss) structures:

Mechanical Engineering Department, GTU

Mechanical Engineering Department, GTU

Stiffness Matrix - Direct Method

Assuming that the displacement u is varying linearly along the axis of

: elongation, F: force in bar

Mechanical Engineering Department, GTU

Stiffness Matrix - A Formal Approach

B: Element strain-displacement matrix

Mechanical Engineering Department, GTU

Stiffness Matrix - A Formal Approach

Stress can be written as:

The work done by the two nodal forces is:

For conservative system, we state that:

Mechanical Engineering Department, GTU

Stiffness Matrix - A Formal Approach

k is the element stiffness matrix. The above expression is a general

which is the same as we derived using the direct method.

Mechanical Engineering Department, GTU

Degree of Freedom (dof)

Number of components of the displacement vector at a node.

Mechanical Engineering Department, GTU

Mechanical Engineering Department, GTU

a) Element Stiffness Matrices :

Mechanical Engineering Department, GTU

Construction of global stiffness matrix :

Mechanical Engineering Department, GTU

Mechanical Engineering Department, GTU

For the data given;

Mechanical Engineering Department, GTU

d) Stresses in the elements:

Mechanical Engineering Department, GTU

In this case, the calculated stresses in elements 1 and 2 are exact

Mechanical Engineering Department, GTU

Solution procedure with matlab

As linear bar element has 2 degrees of freedom (1 at each

At this step boundary conditions are applied manually to the

once the unknown

Solution procedure with matlab

Mechanical Engineering Department, GTU

Matlab functions used

Mechanical Engineering Department, GTU

Matlab functions used

Mechanical Engineering Department, GTU

Matlab functions used

Mechanical Engineering Department, GTU

Matlab functions used

Mechanical Engineering Department, GTU

Mechanical Engineering Department, GTU

Step 2-Copying relevant files and starting Matlab

Mechanical Engineering Department, GTU

Start solving the problem in Command Window:

Mechanical Engineering Department, GTU

Step 4-Assembling the global stiffness matrix:

Mechanical Engineering Department, GTU