University of Florida/Eml4507/Team7 Report2

Problem 1

Find

Verify the equilibrium of node 2; Explain the assembly process; Using CALFEM verify the results of the 2-D truss system with 2 inclined members ^[1]

Given

Figure 1.1: Free body diagram of node 2 of the 2-D truss system with 2 inclined members

$P=7;f_{3}^{(1)}=4.4378f_{4}^{(1)}=2.5622;f_{1}^{(2)}=-4.4378;f_{2}^{(2)}2=4.4578$

Solution

Equilibrium of Node 2

${\begin{aligned}&\Sigma F_{x}=0\\\end{aligned}}$

(1.1)

${\begin{aligned}&\left(-f_{3}^{(1)}\right)+0+\left(-f_{1}^{(2)}\right)=0\\\end{aligned}}$

(1.2)

${\begin{aligned}&-4.4378+0+4.4378=0\\\end{aligned}}$

${\begin{aligned}&\Sigma F_{y}=0\\\end{aligned}}$

(1.3)

${\begin{aligned}&\left(-f_{4}^{(1)}\right)+P+\left(-f_{2}^{(2)}\right)=0\\\end{aligned}}$

(1.4)

${\begin{aligned}&-2.5622+7-4.4578=0\\\end{aligned}}$

Explanation for the Assembly Process

To assemble the global stiffness matrix for the structure; sum the force-displacement matrices in global coordinates:

{\begin{Bmatrix}F_{1}\\F_{2}\\F_{3}\\F_{4}\\F_{5}\\F_{6}\\\end{Bmatrix}}={\begin{bmatrix}f_{1}^{(1)}\\f_{2}^{(1)}\\f_{3}^{(1)}\\f_{4}^{(1)}\\f_{5}^{(1)}\\f_{6}^{(1)}\\\end{bmatrix}}+{\begin{bmatrix}f_{1}^{(2)}\\f_{2}^{(2)}\\f_{3}^{(2)}\\f_{4}^{(2)}\\f_{5}^{(2)}\\f_{6}^{(2)}\\\end{bmatrix}}

(1.5)

Substituting for the force displacement matrix for each element:

{\begin{Bmatrix}F_{1}\\F_{2}\\F_{3}\\F_{4}\\F_{5}\\F_{6}\\\end{Bmatrix}}={\begin{bmatrix}k_{11}^{(1)}&k_{12}^{(1)}&k_{13}^{(1)}&k_{14}^{(1)}&k_{15}^{(1)}&k_{16}^{(1)}\\k_{21}^{(1)}&k_{22}^{(1)}&k_{23}^{(1)}&k_{24}^{(1)}&k_{25}^{(1)}&k_{26}^{(1)}\\k_{31}^{(1)}&k_{32}^{(1)}&k_{33}^{(1)}&k_{34}^{(1)}&k_{35}^{(1)}&k_{36}^{(1)}\\k_{41}^{(1)}&k_{42}^{(1)}&k_{43}^{(1)}&k_{44}^{(1)}&k_{45}^{(1)}&k_{46}^{(1)}\\k_{51}^{(1)}&k_{52}^{(1)}&k_{53}^{(1)}&k_{54}^{(1)}&k_{55}^{(1)}&k_{56}^{(1)}\\k_{61}^{(1)}&k_{62}^{(1)}&k_{63}^{(1)}&k_{64}^{(1)}&k_{65}^{(1)}&k_{66}^{(1)}\\\end{bmatrix}}{\begin{Bmatrix}d_{1}\\d_{2}\\d_{3}\\d_{4}\\d_{5}\\d_{6}\\\end{Bmatrix}}+{\begin{bmatrix}k_{11}^{(2)}&k_{12}^{(2)}&k_{13}^{(2)}&k_{14}^{(2)}&k_{15}^{(2)}&k_{16}^{(2)}\\k_{21}^{(2)}&k_{22}^{(2)}&k_{23}^{(2)}&k_{24}^{(2)}&k_{25}^{(2)}&k_{26}^{(2)}\\k_{31}^{(2)}&k_{32}^{(2)}&k_{33}^{(2)}&k_{34}^{(2)}&k_{35}^{(2)}&k_{36}^{(2)}\\k_{41}^{(2)}&k_{42}^{(2)}&k_{43}^{(2)}&k_{44}^{(2)}&k_{45}^{(2)}&k_{46}^{(2)}\\k_{51}^{(2)}&k_{52}^{(2)}&k_{53}^{(2)}&k_{54}^{(2)}&k_{55}^{(2)}&k_{56}^{(2)}\\k_{61}^{(2)}&k_{62}^{(2)}&k_{63}^{(2)}&k_{64}^{(2)}&k_{65}^{(2)}&k_{66}^{(2)}\\\end{bmatrix}}{\begin{Bmatrix}d_{1}\\d_{2}\\d_{3}\\d_{4}\\d_{5}\\d_{6}\\\end{Bmatrix}}

(1.6)

Substituting zero for where there is no stiffness for each element and changing to local coordinates:

{\begin{Bmatrix}F_{1}\\F_{2}\\F_{3}\\F_{4}\\F_{5}\\F_{6}\\\end{Bmatrix}}={\begin{bmatrix}k_{11}^{(1)}&k_{12}^{(1)}&k_{13}^{(1)}&k_{14}^{(1)}&0&0\\k_{21}^{(1)}&k_{22}^{(1)}&k_{23}^{(1)}&k_{24}^{(1)}&0&0\\k_{31}^{(1)}&k_{32}^{(1)}&k_{33}^{(1)}&k_{34}^{(1)}&0&0\\k_{41}^{(1)}&k_{42}^{(1)}&k_{43}^{(1)}&k_{44}^{(1)}&0&0\\0&0&0&0&0&0\\0&0&0&0&0&0\\\end{bmatrix}}{\begin{Bmatrix}d_{1}\\d_{2}\\d_{3}\\d_{4}\\d_{5}\\d_{6}\\\end{Bmatrix}}+{\begin{bmatrix}0&0&0&0&0&0\\0&0&0&0&0&0\\0&0&k_{11}^{(2)}&k_{12}^{(2)}&k_{13}^{(2)}&k_{14}^{(2)}\\0&0&k_{21}^{(2)}&k_{22}^{(2)}&k_{23}^{(2)}&k_{24}^{(2)}\\0&0&k_{31}^{(2)}&k_{32}^{(2)}&k_{33}^{(2)}&k_{34}^{(2)}\\0&0&k_{41}^{(2)}&k_{42}^{(2)}&k_{44}^{(2)}&k_{44}^{(2)}\\\end{bmatrix}}{\begin{Bmatrix}d_{1}\\d_{2}\\d_{3}\\d_{4}\\d_{5}\\d_{6}\\\end{Bmatrix}}

(1.7)

Summing the matrices:

{\begin{Bmatrix}F_{1}\\F_{2}\\F_{3}\\F_{4}\\F_{5}\\F_{6}\\\end{Bmatrix}}={\begin{bmatrix}k_{11}^{(1)}&k_{12}^{(1)}&k_{13}^{(1)}&k_{14}^{(1)}&0&0\\k_{21}^{(1)}&k_{22}^{(1)}&k_{23}^{(1)}&k_{24}^{(1)}&0&0\\k_{31}^{(1)}&k_{32}^{(1)}&\left(k_{33}^{(1)}+k_{11}^{(2)}\right)&\left(k_{34}^{(1)}+k_{12}^{(2)}\right)&k_{13}^{(2)}&k_{14}^{(2)}\\k_{41}^{(1)}&k_{42}^{(1)}&\left(k_{43}^{(1)}+k_{21}^{(2)}\right)&\left(k_{44}^{(1)}+k_{22}^{(2)}\right)&k_{23}^{(2)}&k_{24}^{(2)}\\0&0&k_{31}^{(2)}&k_{32}^{(2)}&k_{33}^{(2)}&k_{34}^{(2)}\\0&0&k_{41}^{(2)}&k_{42}^{(2)}&k_{43}^{(2)}&k_{44}^{(2)}\\\end{bmatrix}}{\begin{Bmatrix}d_{1}\\d_{2}\\d_{3}\\d_{4}\\d_{5}\\d_{6}\\\end{Bmatrix}}

(1.8)

CALFEM for 2-D Truss System

MATLAB script file:

 
Edof=[1 1 2 3 4;
      2 3 4 5 6;];
bc=[1 0;2 0;5 0;6 0;]
f=zeros(6,1);
f(4)=7;
K=zeros(6);
x1=-2*sqrt(3); y1=-2;
x2=0; y2=0;
x3=sqrt(2); y3=-sqrt(2);
A1=1; A2=2;
E1=3; E2=5;
ex1=[x1 x2]; ey1=[y1 y2]; ep1=[E1 A1];
ex2=[x2 x3]; ey2=[y2 y3]; ep2=[E2 A2];

Ke1=bar2e(ex1,ey1,ep1)
Ke2=bar2e(ex2,ey2,ep2)

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

a=solveq(K,f,bc);

ed1=extract(Edof(1,:),a)
ed2=extract(Edof(2,:),a)

f1=Ke1*ed1'
f2=Ke2*ed2'

MATLAB output:

 
>> r2p1

bc =

     1     0
     2     0
     5     0
     6     0


Ke1 =

    0.5625    0.3248   -0.5625   -0.3248
    0.3248    0.1875   -0.3248   -0.1875
   -0.5625   -0.3248    0.5625    0.3248
   -0.3248   -0.1875    0.3248    0.1875


Ke2 =

    2.5000   -2.5000   -2.5000    2.5000
   -2.5000    2.5000    2.5000   -2.5000
   -2.5000    2.5000    2.5000   -2.5000
    2.5000   -2.5000   -2.5000    2.5000


K =

    0.5625    0.3248   -0.5625   -0.3248         0         0
    0.3248    0.1875   -0.3248   -0.1875         0         0
   -0.5625   -0.3248    3.0625   -2.1752   -2.5000    2.5000
   -0.3248   -0.1875   -2.1752    2.6875    2.5000   -2.5000
         0         0   -2.5000    2.5000    2.5000   -2.5000
         0         0    2.5000   -2.5000   -2.5000    2.5000


ed1 =

         0         0    4.3520    6.1271


ed2 =

    4.3520    6.1271         0         0


f1 =

   -4.4378
   -2.5622
    4.4378
    2.5622


f2 =

   -4.4378
    4.4378
    4.4378
   -4.4378

Problem 2

Find

Homogenous L2-ODE-CC, damping ratio zeta, plot of the solution, log decrement delta, quality factor Q and loss factor eta.

Given

Consider the spring-mass damper system on p.53-2.

 $\lambda _{1}=-0.5+i2,\lambda _{2}=-0.5-i2$

Initial Conditions:

 $y(0)=1,y'(0)=0$

Solution

Case: Underdamped

$\Delta =a^{2}-4b<0$
$\alpha =-0.5,\beta =2$
$(\lambda -(-0.5+2i))(\lambda -(-0.5-2i))=0$
$\lambda ^{2}+1\lambda +2.25=0$

$a=1,b=2.25$

$\omega ^{2}=b$

(2.1)

$2\omega \zeta =a$

(2.2)

Solving equation (2.1) and (2.2)
$\omega =1.5$
$\zeta =0.333$

Therefore L2-ODE-CC: $y_{h}''+y_{h}'+2.25y_{h}=0$

$y=C_{1}e^{(\alpha x)}cos(\beta x)+C_{2}e^{(\alpha x)}sin(\beta x)$

$y=C_{1}e^{(-0.5x)}cos(2x)+C_{2}e^{(-0.5x)}sin(2x)$

(2.3)

$y'=-0.5C_{1}e^{(-0.5x)}cos(2x)-2C_{1}e^{(-0.5x)}sin(2x)-0.5C_{2}e^{(-0.5x)}sin(2x)+2C_{2}e^{(-0.5x)}cos(2x)$

(2.4)

Using the two initial conditions and solving equation (2.3) and (2.4) the coefficients are as follows: $C_{1}=1,C_{2}=.25$

Therefore,
$y=e^{(-0.5x)}cos(2x)+0.25e^{(-0.5x)}sin(2x)$

The next step in analyzing the model is comparing and evaluating Log Decrements. This can be done numerically and by a formula. The log decrement in this case was done both ways and then compared to verify the value.
Numerical method: $\delta ^{(n)}=log[(y_{h}^{(n)})/(y_{h}^{(n+1)})]$

The log decrement was calculated by finding the average value: $\delta =[1/3](\delta ^{(1)}+\delta ^{(2)}+\delta ^{(3)})$

Verifying the log decrement value by using the formula: $\delta =(2\pi \zeta )/{\sqrt {1-\zeta ^{2}}}$

Comparing the two values: $\delta =2.21$

Finding the quality factor : $Q=\delta /\pi$
$Q=0.703$

Finding the loss factor: $\eta =1/Q$
$\eta =1.42$

Problem 3

Find

Determine the member force and axial stress in each member of the truss shown in the figure using one of the FE analysis programs in the Appendix. Assume that Young’s modulus is 104 psi and all cross-sections are circular with a diameter of 2 in. Compare the results with the exact solutions that are obtained from the free-body diagram.^[2]

Write a matlab program to assemble the global stiffness matrix, compute the unknown disp dofs, the reaction forces, the member forces; compare the results with exact hand calculation (Statics, Mechanics of Materials).

Run CALFEM to verify the results, but now do this problem in a different way (as in commercial FE codes).

First, establish 2 arrays: (1) the global node coordinate array “coord”, and (2) the elem connectivity array “conn”. next, write a matlab code to loop over the number of element to call the function “bar2e” of CALFEM to compute the element stiffness matrices. with this method, you let your matlab code construct the arrays ex1, ey1, ex2, ey2, etc. for you, using the info from the arrays “coord” and “conn”.

Given

Solution

By Hand Calculations

$\sum F_{x}=0=400lb+T_{1x}$

(3.1)

$T_{1x}=-400lb$

$\sum M_{1}=0=-1200lb(24')-400lb(9')+T_{4y}(36')$

(3.2)

$T_{4y}=900lb$

$\sum F_{y}=0=-1200lb+900lb+T_{1y}$

(3.3)

$T_{1y}=300lb$

$\sum M_{6}=0=-300lb(12')-400lb(9')+F_{8}(9')$

(3.4)

$F_{8}=800lb$

$\sum M_{3}=0=-300lb(24')+F_{2}(9')$

(3.5)

$F_{8}=800lb$

$\sum F_{y}=0=300lb-{\frac {3}{5}}F_{5}$

(3.6)

$F_{5}=500lb$

By MATLAB

MATLAB script file:

 
%X= [element number; EA/L; l; m];
X = [1:9; 174.44 218.06 174.44 290.74 174.44 290.74 218.06 218.06 218.06; 0.80 -1.00 -0.80 0.00 -0.80 0.00 1.00 1.00 1.00; 
     0.60 0.00 0.60 1.00 0.60 1.00 0.00 0.00 0.00]
i=1;
for i=1:9
str=fprintf('For element %d\n',i);
k=X(2,i).*[X(3,i)^2 X(3,i)*X(4,i) -(X(3,i)^2) -X(3,i)*X(4,i); X(3,i)*X(4,i) X(4,i)^2 -X(3,i)*X(4,i) -(X(4,i)^2); 
  -(X(3,i)^2)-X(3,i)*X(4,i) X(3,i)^2 X(3,i)*X(4,i); -X(3,i)*X(4,i) -(X(4,i)^2) X(3,i)*X(4,i) X(4,i)^2]
i=i+1;
end

MATLAB output:

>> report2

X =

   1.0000    2.0000    3.0000    4.0000    5.0000    6.0000    7.0000    8.0000    9.0000
 174.4400  218.0600  174.4400  290.7400  174.4400  290.7400  218.0600  218.0600  218.0600
   0.8000   -1.0000   -0.8000         0   -0.8000         0    1.0000    1.0000    1.0000
   0.6000         0    0.6000    1.0000    0.6000    1.0000         0         0         0

For element 1

k =

 111.6416   83.7312 -111.6416  -83.7312
  83.7312   62.7984  -83.7312  -62.7984
-111.6416  -83.7312  111.6416   83.7312
 -83.7312  -62.7984   83.7312   62.7984

For element 2

k =

 218.0600         0 -218.0600         0
        0         0         0         0
-218.0600         0  218.0600         0
        0         0         0         0

For element 3

k =

 111.6416  -83.7312 -111.6416   83.7312
 -83.7312   62.7984   83.7312  -62.7984
-111.6416   83.7312  111.6416  -83.7312
  83.7312  -62.7984  -83.7312   62.7984

For element 4

k =

        0         0         0         0
        0  290.7400         0 -290.7400
        0         0         0         0
        0 -290.7400         0  290.7400

For element 5

k =

 111.6416  -83.7312 -111.6416   83.7312
 -83.7312   62.7984   83.7312  -62.7984
-111.6416   83.7312  111.6416  -83.7312
  83.7312  -62.7984  -83.7312   62.7984

For element 6

k =

        0         0         0         0
        0  290.7400         0 -290.7400
        0         0         0         0
        0 -290.7400         0  290.7400

For element 7

k =

 218.0600         0 -218.0600         0
        0         0         0         0
-218.0600         0  218.0600         0
        0         0         0         0

For element 8

k =

 218.0600         0 -218.0600         0
        0         0         0         0
-218.0600         0  218.0600         0
        0         0         0         0

For element 9

k =

 218.0600         0 -218.0600         0
        0         0         0         0
-218.0600         0  218.0600         0
        0         0         0         0

Problem 4

Part 1

Find

For each case, derive the corresponding standard L2-ODE-CC, find the overall homogeneous solution, and plot the solution in a separate figure; then plot all 3 solutions in the same figure for comparison.

Given

Two Roots: $\lambda _{1}=-0.5,\lambda _{2}=-1.5$

Initial conditions y(0)=1, y'(0)=0

Solution

Since \lambda_1 and \lambda_2 are unequal and negative, the system is overdamped. Currently we have:

$(\lambda +0.5)(\lambda +1.5)=0$

Standard form

$(\lambda ^{2}+2*\lambda +0.75)=0$

The standard homogeneous L2-ODE-CC, or the equation of motion, is therefore:

$y_{h}''+2y_{h}'+0.75y_{h}=0$

The general solution is given by:

$y_{h}(t)=c_{1}e^{(-0.5t)}+c_{2}e^{(-1.5t)}$

Using the given initial conditions, we solve for the constants:

$y^{'}(0)=0=-0.5c_{1}-1.5c_{2}$

$c_{1}={\frac {3}{4}},c_{2}=-{\frac {1}{4}}$

Plugging these in to find the overall homogeneous solution:

$y_{p}(t)={\frac {3}{4}}e^{(-0.5t)}-{\frac {1}{4}}e^{(-1.5t)}$

The graph of the solution is:

Part 3

Find

Derive the corresponding standard L2-ODE-CC, find the overall homogeneous solution, and plot the solution in a separate figure; then plot all 3 solutions in the same figure for comparison.

Given

Two roots: $\lambda _{1}=-0.5+i2,\lambda _{2}=-1.5-i2$

Initial conditions:

$y(0)=1,y'(0)=v_{0}=0$

Solution

From the characteristic equation:

$(\lambda -\lambda _{1})(\lambda -\lambda _{2})=0$

$(\lambda +0.5-i2)(\lambda +0.5+i2)=0$

Multiplying out to the standard form:

$\lambda ^{2}+\lambda +4.25=0$

$\lambda ^{2}+a\lambda +b=0$

Where:

$a=c/m=1andb=k/m=4.25$

The standard homogeneous L2-ODE-CC, or the equation of motion, is therefore:

$y_{h}''+y_{h}'+4.25y_{h}=0$

The general solution is given by:

$y_{h}(t)=Ae^{-at/2}cos(\omega _{d}-\phi )$

The variables for this equation are found to be:

$\omega _{n}={\sqrt {b}}=2.06<math>rad/s<math>\zeta =a/(2\omega _{n})=0.243$

$\omega _{d}=\omega _{n}{\sqrt {1-\zeta ^{2}}}=2.00$ rad/s

$A={\sqrt {((v_{0}+\zeta \omega _{n}y_{0})^{2}+(y_{0}\omega _{d})^{2})/(\omega _{d})^{2}}}=1.03$

$\phi =tan^{-1}(y_{0}\omega _{d})/(v_{0}+\zeta \omega _{n}y_{0}))=75.9$

Therefore:

$y_{h}(t)=1.03e^{-0.5t}cos(2t-75.9)$

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

Problem 5

Find

The general procedure in linear finite element calculations.

Generate element matrices, assemble element matrices into the Global System of equations. Solve the global system of equations and evaluate element forces.

Given

Figure 5.1: The figure above consists of a system of three linear elastic springs, and the corresponding finite element model. The system of springs is fixed at the ends and loaded by a single load $F$ ^[4]

Solution

Using Matlab, we first define the topology matrix by Edof containing element numbers and global degrees of freedom

>> Edof=[1 1 2; 2 2 3; 3 2 3]
    1     1     2
    2     2     3
    3     2     3

Now the global stiffness matrix K $(3x3)$ of zeros.

>> K=zeros(3,3)
    0     0     0
    0     0     0
    0     0     0

The load vector $f(3x1)$ with load $F=100$ and load it in position 2.

>> f=zeros(3,1); f(2)=100
    0
  100
    0

Stiffness matrices are generated by using the function spring1e
The function is given by:

function [Ke]=spring1e(ep);
k = ep;  
Ke = [ k -k;
     -k  k];

Now we can set the element property as ep for the springs with their respective spring constants $k$ and $2k$ Where $k=1500$

Commands in matlab are as follows:

>> k=1500; ep1=k; ep2=2*k;
>> Ke1=spring1e(ep1)
Ke1 =
       1500       -1500
      -1500        1500
>> Ke2=spring1e(ep2)
Ke2 =
       3000       -3000
      -3000        3000

With the assem function we can assemble stiffness matrices into the global stiffness matrix K. Which leads too:

assem function:

function [K,f]=assem(edof,K,Ke,f,fe)
   [nie,n]=size(edof);
   t=edof(:,2:n);
   for i = 1:nie
     K(t(i,:),t(i,:)) = K(t(i,:),t(i,:))+Ke;
     if nargin==5
        f(t(i,:))=f(t(i,:))+fe;
     end
   end

Stiffness matrix K:

>> K=assem(Edof(1,:),K,Ke2)
K =
       3000       -3000           0
      -3000        3000           0
          0           0           0
>> K=assem(Edof(2,:),K,Ke1)
K =
       3000       -3000           0
      -3000        4500       -1500
          0       -1500        1500
>> K=assem(Edof(3,:),K,Ke2)
K =
       3000       -3000           0
      -3000        7500       -4500
          0       -4500        4500

The Global system of equations is solved with boundary condtions given in bc as well as with function solveq.

Function solveq is as follows:

function [d,Q]=solveq(K,f,bc)
 if nargin==2 ; 
    d=K\f ; 
 elseif nargin==3;
    [nd,nd]=size(K);
    fdof=[1:nd]';
    d=zeros(size(fdof));
    Q=zeros(size(fdof));
    pdof=bc(:,1);
    dp=bc(:,2);
    fdof(pdof)=[];
    s=K(fdof,fdof)\(f(fdof)-K(fdof,pdof)*dp);
    d(pdof)=dp;
    d(fdof)=s;
 end  
    Q=K*d-f;

With Matlab commands it results in:

>> bc=[1 0; 3 0];
>> [a,r]=solveq(K,f,bc)
a =
        0
   0.0133
        0
r =
 -40.0000
        0
 -60.0000

Element forces are now computed from the element displacements. These are obtained from global displacements a using the function extract.

Function extract:

function [ed]=extract(edof,a)
   [nie,n]=size(edof);
   t=edof(:,2:n);
   for i = 1:nie
       ed(i,1:(n-1))=a(t(i,:))';
   end

Now the corresponding Matlab commands are as follows:

>> ed1=extract(Edof(1,:),a)
ed1 =
        0    0.0133
>> ed2=extract(Edof(2,:),a)
ed2 =
   0.0133         0
>> ed3=extract(Edof(3,:),a)
ed3 =
   0.0133         0

The last step is to evaluate the spring forces by using the function spring1s.

We have function spring1s as written below:

function [es]=spring1s(ep,ed)
 k = ep;  
 u=ed;
 es = k*(u(2)-u(1));

We can now finally evaluate the spring forces. Using Matlab the forces are as follows:

>> es1=spring1s(ep2,ed1)
es1 =
   40
>> es2=spring1s(ep1,ed2)
es2 =
  -20
>> es3=spring1s(ep2,ed3)
es3 =
  -40

Note: Instructions and commands were all followed from CALFEM 3.4 Manual.

Problem 6

Find

Find the displacements of the three bodies and the forces ( tensile/compressive) in the springs. What are the reactions at the walls? Assume only horizontal translation.

Given

Solution

The nodes and elements are labeled in the problem. Create a connectivity table.

Construct Elemental equilibrium equations for all 6 elements.

Here is a sample equilibrium equation for element 3:

${\binom {f_{2}^{(3)}}{f_{3}^{(3)}}}=k_{3}{\begin{bmatrix}1&-1\\-1&1\\\end{bmatrix}}{\begin{bmatrix}u_{2}\\u_{3}\\\end{bmatrix}}$

(6.1)

The notation for : The superscript parenthesis signifies the element number, while the subscript signifies the associated node number.

With element forces calculated, a node equilibrium can be established. Start with a free-body diagram of each node, and then create the five nodal equations:

Make sure all forces are positive and modeled in tension.

Combine your elemental equilibrium equations with your nodal equilibrium equations for the global stiffness correlation:

$[F_{s}]=[Q_{s}][K_{s}]$

(6.2)

${\begin{bmatrix}F_{1}=R_{1}\\F_{2}=0\\F_{3}=1000\\F_{4}=0\\F_{5}=R_{5}\end{bmatrix}}={\begin{bmatrix}k_{1}+k_{4}&-k_{1}&-k_{4}&0&0\\-k_{1}&k_{1}+k_{2}+k_{3}&-k_{3}&-k_{2}&0\\-k_{4}&-k_{3}&k_{3}+k_{4}+k_{5}&-k_{5}&0\\0&-k_{4}&-k_{4}&k_{2}+k_{5}+k_{6}&0\\0&0&0&-k_{6}&k_{6}\end{bmatrix}}{\begin{bmatrix}u_{1}\\u_{2}\\u_{3}\\u_{4}\\u_{5}\\\end{bmatrix}}$

(6.3)

Solving via MATLAB

%Spring constants - in N/mm
k1 = 500;
k2 = 400;
k3 = 600;
k4 = 200;
k5 = 400;
k6 = 300;

%Declaration of Global Displacement Variables
u1 = 0; %due to being at reactions
u2 = 1; %placeholder
u3 = 1; %placeholder
u4 = 1; %placeholder
u5 = 0; %due to being at reactions

%Combining the Nodal equilibrium and the Element Equilibrium, we can find
%the displacements

kglobal =[k1+k4   -k1     -k4     0     0;
           -k1  k1+k2+k3  -k3    -k2    0;
           -k4    -k3   k3+k4+k5 -k5    0;
            0     -k2     -k5 k2+k5+k6 -k6;
            0      0       0    -k6     k6;]

%Rows 1 and 5 can be excluded for now, because the displacements at node 1
%and 5 are zero. This makes the corresponding rows and columns irrelevant

kglobalsimplified = [k1+k2+k3  -k3    -k2;
                     -k3     k3+k4+k5 -k5;
                     -k2       -k5  k2+k5+k6;]
                  
%kglobalsimplified*displacements = External Forces

Fsimplified = [0; %no external force
               1000; %given
               0;]; %no external force
           
Displacementsimplified = kglobalsimplified\Fsimplified %u2, u3, u4

u2 = Displacementsimplified(1)
u3 = Displacementsimplified(2)
u4 = Displacementsimplified(3)

klocal = [1 -1;
          -1 1;]; %Local Equilibrium equations are multiplied by this
      
%Element Equilibrium Equations
k1ef = k1*klocal*[u1;u2]; %k1ef means K1 Element Forces
k2ef = k2*klocal*[u2;u4];
k3ef = k3*klocal*[u2;u3];
k4ef = k4*klocal*[u1;u3];
k5ef = k5*klocal*[u3;u4];
k6ef = k6*klocal*[u4;u5];

%The above matrices are 2x1s
%We need to access the specific values within each matrix
%So we'lll reassign the matrices to new variables

f11 = k1ef(1); %k11 = k(node 1)(element 1)
f21 = k1ef(2);

f22 = k2ef(1); 
f42 = k2ef(2);

f23 = k3ef(1); 
f33 = k3ef(2);

f14 = k4ef(1); 
f34 = k4ef(2);

f35 = k5ef(1); 
f45 = k5ef(2);

f46 = k6ef(1); 
f56 = k6ef(2);

%With our element force variables, we can create the nodal equilibirum eqns
F3 = 1000;

R1 = f11 + f14 %Reaction Force
F2 = f21 + f23 + f22; %unknown force
F3 = f33 + f34 + f35; %Given force
F4 = f42 + f45 + f46; %unknown force
R5 = f56 %Reaction Force

Resulting in the answers (in units mm for the displacements and N for the Reaction):

>> run r26.m

u2 =

    0.8542


u3 =

    1.5521


u4 =

    0.8750


R1 =

 -737.5000


R5 =

 -262.5000

Solving via CALFEM

ep1 = 500; %spring constants
ep2 = 400;
ep3 = 600;
ep4 = 200;
ep5 = 400;
ep6 = 300;

%create element spring stiffness matrix
k1 = spring1e(ep1);
k2 = spring1e(ep2);
k3 = spring1e(ep3);
k4 = spring1e(ep4);
k5 = spring1e(ep5);
k6 = spring1e(ep6);

%define degrees of freedom for each element. All elements have 2 dof.
Edof = [1 1 2; % el1 is element 1
        2 2 4;
        3 2 3;
        4 1 3;
        5 3 4;
        6 4 5;];
fe = zeros(5,1); %declare f, and then populate f with extractions
fe(1) = spring1s(ep1,Edof(1,:));
fe(2) = spring1s(ep2,Edof(2,:));
fe(3) = spring1s(ep3,Edof(3,:));
fe(4) = spring1s(ep4,Edof(4,:));
fe(5) = spring1s(ep5,Edof(5,:));

%Assembling structural stiffness matrix
K = zeros(5);
[K,f] = assem(Edof(1,:),K,k1,f,fe);
[K,f] = assem(Edof(2,:),K,k2,f,fe);%adding element 2 to the structural matrix
[K,f] = assem(Edof(3,:),K,k3,f,fe);%adding element 3 to the structural matrix
[K,f] = assem(Edof(4,:),K,k4,f,fe);
[K,f] = assem(Edof(5,:),K,k5,f,fe);
[K,f] = assem(Edof(6,:),K,k6,f,fe);

a = solveq(K,f); %solves K*a = f

Ed = extract(Edof,a) %extracts element displacements from global 
                      %solution vector a
                      
Ed(3)
Ed(4)
Ed(5)

f(1)
f(5)

Resulting in the answers (in units mm for the displacements and N for the Reaction):

>> run calfem26.m

Ed(3) =

    0.8542


Ed(4) =

    1.5521


Ed(5) =

    0.8750


f(1) =

 -737.5000


f(5) =

 -262.5000

References

↑ https://docs.google.com/document/d/1lQLCUXKTYtMs4BqzLExS5jqLVcpAIP1FW4sd2sl1cqA/edit#
↑ KS 2008 p.103 sec.2.7 pb.21
↑ KS 2008 p.103 sec.2.7 pb.21
↑ http://www.solid.lth.se/fileadmin/hallfasthetslara/utbildning/kurser/FHL064_FEM/calfem34.pdf
↑ KS 2008 p.98 sec.2.7 pb.21
↑ http://upload.wikimedia.org/wikiversity/en/0/04/Eml4500.f08.1.djvu

Contributing Team Members

Problem Assignments
Problem #	Solved & Typed by	Reviewed by
1	Matthew Gidel	All
2	Zeyn Kermani	All
3	Kristin Howe	All
4	Spencer Herran,Brandon Wright	All
5	Kevin Li	All
6	Joshua Plicque	All

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

[1] ttps://docs.google.com/document/d/1lQLCUXKTYtMs4BqzLExS5jqLVcpAIP1FW4sd2sl1cqA/edit#

[2] KS 2008 p.103 sec.2.7 pb.21

[3] KS 2008 p.103 sec.2.7 pb.21

[4] ttp://www.solid.lth.se/fileadmin/hallfasthetslara/utbildning/kurser/FHL064_FEM/calfem34.pdf

[5] KS 2008 p.98 sec.2.7 pb.21

[6] ttp://upload.wikimedia.org/wikiversity/en/0/04/Eml4500.f08.1.djvu

[1]

[2]

[3]

[4]

[5]

[6]