User:Eml4507.s13.team3.buckles

= Problem 4: Project 2.1 p.94 sec.2.6 KS 2008 using Matlab and CALFEM = On my honor, I have neither given nor received unauthorized aid in doing this assignment.

Given
Space truss as shown that consists of 25 elements. Each element has a circular cross-section that initially has a diameter $$d=2\;in$$. At the nodes 1 and 2 a constant force $$F=60,000\;lb$$ is applied in the y-direction. The nodes 7, 8, 9 and 10 are fixed. The elements are made of steel with the following properties, Young's modulus $$E=3\times 10^7\; psi$$, Poisson's ratio $$\nu=0.3$$, yield stress $$\sigma_y=37,000\; psi$$ and mass density $$\rho=7.3\times10^-4\;\frac{lb*sec^2}{in^4}$$. A safety factor $$N=1.5$$ is used in the design of the truss. The diameter of the truss elements should be between 0.1 in and 2.5 in.

Find
1. Solve the initial truss structure problem using matlab. Provide plots of the deformed truss and a table of the stresses in the beams.

2. Minimize the structural weight by changing the cross-sectional diameter for each truss element. Identify the zero force members and provide a plot of the deformed shape under optimal design. Calculate the weights of the original truss and the optimal truss.

3. Verify the problem with CALFEM.

1. Original Truss Deformation Using Matlab
Enter the number of global nodes. >> G=10;

Enter the x, y, and z coordinates. >> x=[-1.5 1.5 -1.5 1.5 1.5 -1.5 4 4 -4]; >> y=[0 0 1.5 1.5 -1.5 -1.5 4 4 -4 -4]; >> z=[8 8 4 4 4 4 0 0 0 0];

Enter the input force matrix. The first element of each row corresponds to the global degree of freedom; the second element corresponds to the force associated with that degree of freedom. Only the known forces are entered into this matrix. >> F_in=[1 0; 2 60000; 3 0; 4 0; 5 60000; 6 0; 7 0; 8 0; 9 0; 10 0; 11 0; 12 0; 13 0; 14 0; 15 0; 16 0; 17 0; 18 0];

Enter the input displacement matrix. The first element of each row corresponds to the global degree of freedom; the second element corresponds to the displacement associated with that degree of freedom. >> Q_in=[19 0; 20 0; 21 0; 22 0; 23 0; 24 0; 25 0; 26 0; 27 0; 28 0; 29 0; 30 0];

Enter the number of elements. >> E=25;

Enter the global node corresponding to the first local node of each element. The column number corresponds to the element number of the given node. >> l1=[1 1 2 1 2 2 2 1 1 3 4 3 5 3 6 4 5 4 3 5 6 6 3 4 5];

Enter the global node corresponding to the second local node of each element. >> l2=[2 4 3 5 6 4 5 3 6 6 5 4 6 10 7 9 8 7 8 10 9 10 7 8 9];

Enter the Young's Modulus associated with each element. The column number corresponds to the element number of the given node. >> Modulus=3e7*[1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1];

Enter the Area associated with each element. The column number corresponds to the element number of the given node. >> Area=(pi/4)*2^2*[1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1];

Enter the Mass Density associated with each element. The column number corresponds to the element number of the given node. >> Ro=7.3e-4*[1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1]; Run the truss program. >> [F_bar] = truss(G, x, y, z, F_in, Q_in, E, l1, l2, Modulus, Area, Ro)

The internal spring force matrix F_bar F_bar = 1.0e+04 * Columns 1 through 9 0.0000   3.5997    3.5997   -3.5997   -3.5997    5.5981   -5.5981    5.5981   -5.5981      -0.0000   -3.5997   -3.5997    3.5997    3.5997   -5.5981    5.5981   -5.5981    5.5981     Columns 10 through 18 -0.0000  -0.0000   -0.9053    0.9053    1.8114   -1.8114    1.8114   -1.8114    3.4748       0.0000    0.0000    0.9053   -0.9053   -1.8114    1.8114   -1.8114    1.8114   -3.4748     Columns 19 through 25 3.4748  -3.4748   -3.4748   -6.7822    6.7822    6.7822   -6.7822      -3.4748    3.4748    3.4748    6.7822   -6.7822   -6.7822    6.7822

Row 1 gives the member forces with respect to the Global nodes and row 2 gives the member forces with respect to the members.

The axial stresses are given by dividing the axial forces by the area. stress = 1.0e+04 * Columns 1 through 9 -0.0000  -1.1458   -1.1458    1.1458    1.1458   -1.7819    1.7819   -1.7819    1.7819     Columns 10 through 18 0.0000   0.0000    0.2882   -0.2882   -0.5766    0.5766   -0.5766    0.5766   -1.1061     Columns 19 through 25 -1.1061   1.1061    1.1061    2.1588   -2.1588   -2.1588    2.1588



Deformed (red) and undeformed (black) plot of the Truss. The shown deformation is scaled 10 times the actual deformation.

The weight, which is found by summing up the mass matrix and multiplying by the acceleration due to gravity and dividing by 3 is 117.17 lb.

Matlab Code
This code is equivalent to that used in R2.3 except that it uses the given spring constants rather than than the modulus and area to determine the stiffness matrix.

2. Optimal Truss Deformation Using Matlab
In order to minimize the diameter of the elements the stress in each element was used with the yield stress, a factor of safety, and Poisson's ratio to construct a formula for the new diameter of the beams.

New Diameter Calculations Stress in a beam. $$\sigma =\frac{P}{A}$$ Yield stress scaled with a factor of safety. $$\frac{\sigma_{y}}{1.5} =\frac{P}{A}$$ Solving for deformed diameter. $$D=\sqrt{\frac{4}{\pi}\frac{P}{\sigma_{y}/N}}$$ Poisson's ratio. $$\nu = -\frac{d\varepsilon_\mathrm{trans}}{d\varepsilon_\mathrm{axial}} = -\frac{d\varepsilon_\mathrm{y}}{d\varepsilon_\mathrm{x}}= -\frac{d\varepsilon_\mathrm{z}}{d\varepsilon_\mathrm{x}}$$ Poisson's ratio solved for a circular beam. The variable "d" is the new undeformed diameter, "l" is the deformed length and "L" is the undeformed length. $$\nu=-\frac{\frac{D-d}{d}}{\frac{l-L}{L}}$$ Poisson's ratio solved for a the new diameter. $$d=\frac{D}{\frac{-0.3(l-L)}{L}+1}$$ Matlab Code that calculates the new diameters. Table of the stress in each member with the new diameters. Conclusion The zero force members are 1, 10 and 11. The new diameters allow for the an approximate stress of 24666 psi under the load given in the problem statement. From the figure of the deformed truss it can be seen that this truss deforms much more than the original because of the smaller areas. However, this truss weighs only 54.35 pounds while the original weighs 117.17 pounds.

3. CALFEM Verification
First begin by constructing the Edof  matrix. Column 1 corresponds to the element number; columns 2 and 3 correspond to the global node that corresponds to the first and second local nodes respectively. >> Edof=[1 1 2 5 6 2 1 2 3 4           3 1 2 7 8];

Construct an empty stiffness matrix K that has side lengths equal to the number of degrees of freedom times the number of global nodes. For this problem, there are 4 global nodes and motion is constrained to 2 degrees of freedom. >> K=zeros(8,8);

Construct the force matrix F. The forces corresponding to global node 1 are known; the forces corresponding to global nodes 2, 3, and 4 are unknown. >> F=zeros(8,1); >> F(1)=20000/sqrt(2); >> F(2)=20000/sqrt(2);

Construct the displacement matrix Q. Global nodes 2, 3, and 4 and known displacements as shown. >> Q=[3 .02 4 -.01        5 -.03         6 .05         7 0         8 0];

The function bar2e generates the element stiffness matrices for a bar with 2 degrees of freedom. Its parameters are the x and y coordinates of the element and a vector containing the Young's modulus and cross sectional area. This information is given in the problem statement.

The Vector containing the Young's modulus and cross sectional area. >> E=206e9; >> A=1e-4; >> ep=[E A];

The x and y coordinates for each element. >> ex1=[0 sqrt(3)/2]; >> ey1=[0 -.5]; >> ex2=[0 0]; >> ey2=[0 1]; >> ex3=[0 -sqrt(3)/2]; >> ey3=[0 -.5];

The commands: >> Ke1=bar2e(ex1, ey1, ep); >> Ke2=bar2e(ex2, ey2, ep); >> Ke3=bar2e(ex3, ey3, ep);

Yield the following stiffness matrices: >> Ke1 = 1.0e+07 * 1.5450  -0.8920   -1.5450    0.8920          -0.8920    0.5150    0.8920   -0.5150          -1.5450    0.8920    1.5450   -0.8920           0.8920   -0.5150   -0.8920    0.5150      Ke2 = 0          0           0           0                0    20600000           0   -20600000                0           0           0           0                0   -20600000           0    20600000   >> Ke3 = 1.0e+07 * 1.5450   0.8920   -1.5450   -0.8920          0.8920    0.5150   -0.8920   -0.5150         -1.5450   -0.8920    1.5450    0.8920         -0.8920   -0.5150    0.8920    0.5150

The global stiffness matrix can be assembled with the element stiffness matrices using the following command. >> K=assem(Edof(1,:),K,Ke1); >> K=assem(Edof(2,:),K,Ke2); >> K=assem(Edof(3,:),K,Ke3);

This results in the following stiffness matrix. >> K = 1.0e+07 * 3.0900        0         0         0   -1.5450    0.8920   -1.5450   -0.8920               0    3.0900         0   -2.0600    0.8920   -0.5150   -0.8920   -0.5150               0         0         0         0         0         0         0         0               0   -2.0600         0    2.0600         0         0         0         0         -1.5450    0.8920         0         0    1.5450   -0.8920         0         0          0.8920   -0.5150         0         0   -0.8920    0.5150         0         0         -1.5450   -0.8920         0         0         0         0    1.5450    0.8920         -0.8920   -0.5150         0         0         0         0    0.8920    0.5150

This stiffness matrix is in accord with that computed in Matlab in the previous section.

The function solveq takes the stiffness matrix, known force matrix, and displacement matrix and returns the complete force and displacement matrices. >> [q r]=solveq(K,F,Q) q = -0.0290      0.0108       0.0200      -0.0100      -0.0300       0.0500            0            0   r = 1.0e+05 * 0.0000      0.0000            0      -4.2816      -3.6562       2.1109       3.5148       2.0293

The displacement and force matrices are in accord with those computed in Matlab in the previous section.

The function extract is used to evaluate the displacements of the local nodes for each beam, which is needed to determine the force in each beam. >> ed1=extract(Edof(1,:),q); >> ed2=extract(Edof(2,:),q); >> ed3=extract(Edof(3,:),q);

This results in the following displacements. The number following 'ed' corresponds to the element number. Each column corresponds to the respective degree of freedom for the beam. ed1 = -0.0290   0.0108   -0.0300    0.0500 ed2 = -0.0290   0.0108    0.0200   -0.0100 ed3 = -0.0290   0.0108         0         0

The function bar2s can be used to determine the force in each beam. Its input parameters are the x and y coordinates for the beam (ex1 and ey1), the physical properties of the beam (Young's modulus and cross sectional area - ep) and the displacement (ed). >> es1=bar2s(ex1, ey1,ep,ed1) >> es1=bar2s(ex2, ey2,ep,ed2) >> es1=bar2s(ex3, ey3,ep,ed3)

This results in the following spring forces. es1 = -4.2219e+05 es2 = -4.2816e+05 es3 = -4.0586e+05

This is in accord with the forces computed in Matlab in the previous section.