exs_beam1¶

Purpose:

Analysis of a simply supported beam.

Description:

Consider a beam with the length \(9.0\) m. The beam is simply supported and loaded by a point load \(P=10000\) N applied at a point \(3.0\) m from the left support. The corresponding computational model has six degrees of freedom and consists of two beam elements with four degrees of freedom. The beam has Young’s modulus \(E=210\) GPa and moment of inertia \(I=2510 \times 10^{-8}\) m⁴.

Example:

The element topology is defined by the topology matrix:

>>  Edof=[1  1  2  3  4;
          2  3  4  5  6];

The system matrices, i.e. the stiffness matrix K and the load vector f, are defined by:

>> K=zeros(6);   f=zeros(6,1);   f(3)=-10000;

The element property vector ep, the element coordinate vectors ex1 and ex2, and the element stiffness matrices Ke1 and Ke2, are generated:

>> E=210e9;     I=2510e-8;      ep=[E A I];
>> ex1=[0 3];     ex2=[3 9];

>> Ke1=beam1e(ex1,ep)

Ke1 =

   1.0e+06 *

    2.3427    3.5140   -2.3427    3.5140
    3.5140    7.0280   -3.5140    3.5140
   -2.3427   -3.5140    2.3427   -3.5140
    3.5140    3.5140   -3.5140    7.0280

>> Ke2=beam1e(ex2,ep)

Ke2 =

   1.0e+06 *

    0.2928    0.8785   -0.2928    0.8785
    0.8785    3.5140   -0.8785    1.7570
   -0.2928   -0.8785    0.2928   -0.8785
    0.8785    1.7570   -0.8785    3.5140

Based on the topology information, the global stiffness matrix can be generated by assembling the element stiffness matrices:

>> K=assem(Edof(1,:),K,Ke1);
>> K=assem(Edof(2,:),K,Ke2);

Finally, the solution can be calculated by defining the boundary conditions in bc and solving the system of equations. Displacements a and support forces r are computed by the function solveq:

>> bc=[1 0; 5 0];
[a,r]=solveq(K,f,bc)

The section forces es are calculated from element displacements Ed:

>> Ed=extract_ed(Edof,a);

>> [es1,edi1]=beam1s(ex1,ep,Ed(1,:),eq,6)

>> [es2,edi2]=beam1s(ex2,ep,Ed(2,:),eq,11)

Results:

a =                      r =
         0                   1.0e+003 *
   -0.0095
   -0.0228                     6.6667
   -0.0038                    -0.0000
         0                    -0.0000
    0.0076                    -0.0000
                               3.3333
                                    0

es1 =                        edi1 =

  1.0e+004 *                           0
                                 -0.0093
   -0.6667    0.0000             -0.0173
   -0.6667    0.6667             -0.0228
   -0.6667    1.3333
   -0.6667    2.0000

es2 =                        edi2 =

 1.0e+004 *                      -0.0228
                                 -0.0248
    0.3333    2.0000             -0.0236
    0.3333    1.6667             -0.0199
    0.3333    1.3333             -0.0143
    0.3333    1.0000             -0.0075
    0.3333    0.6667             -0.0000
    0.3333    0.3333
    0.3333   -0.0000

A displacement diagram and section force diagrams are displayed using the function plot:

figure(1)
hold on;
plot([0 9],[0 0]);
c=plot([0,0:1:3,3:1:9,9],[0;edi1(:,1);edi2(:,1);0]);
set(c,'LineWidth',[2]);
axis([-1 10 -0.03 0.01]);
title('displacements')

figure(2)
hold on;
plot([0 9],[0 0]);
c=plot([0,0:1:3,3:1:9,9],[0;es1(:,1);es2(:,1);0]);
set(c,'LineWidth',[2]);
axis([-1 10 -8000 5000]);
set(gca, 'YDir','reverse');
title('shear force')

figure(3)
hold on;
plot([0 9],[0 0]);
c=plot([0,0:1:3,3:1:9,9],[0;es1(:,2);es2(:,2);0]);
set(c,'LineWidth',[2]);
axis([-1 10 -5000 25000]);
set(gca, 'YDir','reverse');
title('moment')