User:Eml4500.f08.RAMROD.A/HW4

 One member needed a little extra time to complete his section: the section titled, "The Transfer Function," so it was uploaded by him and then placed in the right position by the team coordinator. Although this comment box was added today (Eml4500.f08.RAMROD.A 04:37, 29 October 2008 (UTC)) the two versions of HW 4 were completed within a few minutes. One version is on time and the other approximately 13 minutes late. [Here is the Comparison]

Connectivity Array "conn"
The following array is the array that should be used in Matlab to represent the global node numbers corresponding to each element of the truss system. The rows of this array represent the element number and the columns represent the local node number.

$$\displaystyle conn = \begin{bmatrix} 1 & 2\\ 2 & 3 \end{bmatrix}$$

Row one of this array corresponds to element one while row two corresponds to element two. Column one represents local node one and column two represents local node two. So to interpret this array we will do an example. Lets say we want to find the global node number for local node one of element one. Since the columns represent the local node number we go to column one and since the rows represent the element number we go to row one. In spot $$\displaystyle (1,1)$$ of the conn array it can be determined that the global node number of element one, local node one is simply one. In general;

$$\displaystyle conn(e,j) = $$ global node number of local node j of element

Location Matrix Master Array "Lmm"
The following array is similar to that of the array above with a minor difference. It is used to determine the global degrees of freedom when the local degrees of freedom and element number are known. The rows represent the element number and instead of representing the local node number, the columns represent the local degrees of freedom. Since the problem being analyzed is a two truss problem, each element has four degrees of freedom so the array has four columns and two rows representing each element.

$$\displaystyle Lmm = \begin{bmatrix} 1 & 2 &3 &4 \\ 3& 4 & 5 & 6 \end{bmatrix}$$

To interpret this array we will do an example. Lets say we want to find the global degree of freedom of element two at local degree of freedom three. With this information we go to row two (element two) and column three (local degree of freedom). At this location we find the global degree of freedom corresponding to element two, local degree of freedom three. In general;

$$\displaystyle Lmm(e,j)=$$ the global degree of freedom number for element $$\displaystyle e$$ corresponding to the $$\displaystyle j^{th}$$ local degree of freedom

The Transfer Function
Goal - We want to find T that transfer the set of local elemental dofs (degrees of freedom), d, to another set of local elemental dofs, d, such that T is invertible.



$$\tilde{\underline{d}}^{(e)}_{4x1} = \tilde{\underline{T}}^{(e)}_{4x4}\;\tilde{\underline{d}}^{(e)}_{4x1}$$

Equation (1): $$\tilde{d}_1^{(e)}=\begin{bmatrix} l^{(e)}&m^{(e)} \end{bmatrix} = \begin{Bmatrix} d_1^{(e)}\\ d_2^{(e)} \end{Bmatrix}$$

$$\tilde{d}_2^{(e)} = \overrightarrow{d}_1^{(e)} \overrightarrow{\tilde{j}}$$

$$\tilde{d}_2^{(e)} = -\sin \theta^{(e)}d_1^{(e)} + \cos \theta^{(e)}d_2^{(e)}$$

Equation (2): $$\tilde{d}_2^{(e)} = \begin{bmatrix} -m^{(e)} & l^{(e)} \end{bmatrix} = \begin{Bmatrix} d_1^{(e)}\\ d_2^{(e)} \end{Bmatrix}$$

Put Equations (1) and (2) together in matrix form:

$$\begin{Bmatrix} \tilde{d}_1^{(e)}\\ \tilde{d}_2^{(e)} \end{Bmatrix} = \begin{bmatrix} l^{(e)} & m^{(e)}\\ -m^{(e)}& l^{(e)} \end{bmatrix} \begin{Bmatrix} d_1^{(e)}\\ d_2^{(e)} \end{Bmatrix}$$

$$\begin{Bmatrix} \tilde{d}_1^{(e)}\\ \tilde{d}_2^{(e)}\\ \tilde{d}_3^{(e)}\\ \tilde{d}_4^{(e)} \end{Bmatrix} = \begin{bmatrix} \underline{R}_{2x2}^{(e)} & \underline{O}_{2x2}\\ \underline{O}_{2x2}& \underline{R}_{2x2}^{(e)} \end{bmatrix} \begin{Bmatrix} d_1^{(e)}\\ d_2^{(e)}\\ d_3^{(e)}\\ d_4^{(e)} \end{Bmatrix}$$



$$\tilde{\underline{f}}^{(e)} = k^{(e)}\begin{vmatrix} 1 &0 &-1  &0 \\ 0 &0  & 0 &0 \\ -1& 0 & 1 & 0\\ 0 & 0 &0  &0 \end{vmatrix} \tilde{d}^{(e)}$$

$$\tilde{\underline{f}}_{4x1}^{(e)} = \tilde{\underline{k}}_{4x4}^{(e)}\;\tilde{\underline{d}}_{4x1}^{(e)}$$

Matrix Diagonals
Note: Consider case where $$\tilde{d}_{4}^{(e)} \ne 0$$

$$\tilde{d}_{1}^{(e)}=\tilde{d}_{2}^{(e)}=\tilde{d}_{3}^{(e)}= 0$$

$$\tilde{\mathbf{f}}_{4x1}^{(e)}=\tilde{\mathbf{k}}_{4x4}^{(e)}\tilde{\mathbf{d}}_{4x1}^{(e)}=\mathbf{0}_{4x1}$$

The $$\mathbf{0}_{4x1}$$ is the 4th column of $$\mathbf{k}^{(e)}$$

Interpretation of the the transverse degrees of freedom

$$\tilde{\mathbf{d}}^{(e)}=\tilde{\mathbf{T}}^{(e)}{d}^{(e)}$$

Similarly: $$\tilde{\mathbf{f}}^{(e)}=\tilde{\mathbf{T}}^{(e)}{f}^{(e)}$$

Also:

$$\tilde{\mathbf{k}}^{(e)}{\tilde{\mathbf{d}}^{(e)}}=\tilde{\mathbf{f}}^{(e)}$$

Therefore:

$$\tilde{\mathbf{k}}^{(e)}\tilde{\mathbf{T}}^{(e)}\mathbf{d}^{(e)}=\tilde{\mathbf{T}}^{(e)}\mathbf{f}^{(e)}$$

If $$\tilde{\mathbf{T}}^{(e)}$$ is invertible, then:

$$ [\tilde{\mathbf{T}}^{(e)-1}{\tilde{\mathbf{k}}^{(e)}}\tilde{\mathbf{T}}^{(e)}]\mathbf{d}^{(e)}=\mathbf{f}^{(e)}$$

$$\tilde{\mathbf{T}}^{(e)} $$ block diagonal matrix: $$\begin{bmatrix} \mathbf{R}_{2x2}^{(e)} & \mathbf{0}_{2x2} \\ \mathbf{0}_{2x2} & \mathbf{R}_{2x2}^{(e)} \end{bmatrix}$$

Consider a general block diagonal matrix: $$ \mathbf{A}=\begin{bmatrix} \mathbf{D}_{1}& \mathbf{0} \\ \mathbf{0} & \mathbf{D}_{s} \end{bmatrix}$$

Question: What is $$\mathbf{A}^{-1}$$?

Here is a simple example

$$ \mathbf{B}=\begin{bmatrix} \mathbf{d}_{11}& \mathbf{0}& \mathbf{0} \\ \mathbf{0} & \mathbf{d}_{22} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{d}_{nn} \end{bmatrix}$$

$$\mathbf{B}=Diag[d_{11}, d_{22}, ..., d_{nn}]$$ (This is the shortened notation for $$ \mathbf{B} $$)

$$ \mathbf{B}^{-1}=Diag[1/d_{11}, 1/d_{22}, ..., 1/d_{nn}]$$

Assumming $$d_{ii}\ne 0$$ for $$ i=1, ..., n$$

For a block diagonal matrix A:

$$\mathbf{A}=Diag[d_{1} ..., d_{s}]$$

$$\mathbf{A}^{-1}=Diag[d_{1}^{-1} ..., d_{s}^{-1}]$$

$$\tilde{\mathbf{T}}^{(e)-1}=Diag[\mathbf{T}^{(e)-1} ..., R^{(e)-1}]$$

We know:

$$\mathbf{R}^{(e)T}=\begin{bmatrix} l^{(e)}& -m^{(e)} \\ m^{(e)} & -l^{(e)} \end{bmatrix}$$

Also

$$\mathbf{R}_{2x2}^{(e)}\mathbf{R}_{2x2}^{(e)T}=\begin{bmatrix} 1& 0 \\ 0 & 1 \end{bmatrix}_{2x2}=\mathbf{I}_{2x2}$$ (This is the Identity matrix)

In conclusion: $$\mathbf{R}^{(e)-1}=\mathbf{R}^{(e)T}$$

So $$\tilde{\mathbf{T}}^{(e)-1}=(Diag[\mathbf{R}^{(e)} ,\mathbf{R}^{(e)}])^{T}$$ where $$Diag[\mathbf{R}^{(e)} ,\mathbf{R}^{(e)}]=\tilde{\mathbf{T}}^{(e)}$$

This shows that: $$\tilde{\mathbf{T}}^{(e)-1}=\tilde{\mathbf{T}}^{(e)T}$$

HW: Eigenvectors for the zero Eigenvalues of K
Using Matlab, the eigenvalues for K were found (see HW 3) and the eigenvectors were plotted for the zero values. The eigenvectors were found to be:

$$E_1 = \begin{Bmatrix} -0.4494 & 0.6754 & 0.1682 & -0.3942 & 0.2812 & -0.2812 \end{Bmatrix}$$

$$E_2 = \begin{Bmatrix} 0.3877&-0.4281&0.3877&-0.4281&0.4163&-0.3995 \end{Bmatrix}$$

$$E_3 = \begin{Bmatrix} 0.5158&0.4825&0.5158&0.4825&-0.0116&-0.0449 \end{Bmatrix}$$

$$E_4 = \begin{Bmatrix} 0.0160&0.0244&0.0160&0.0244&0.7023&0.7107 \end{Bmatrix}$$

The plots are shown below.







These plots show the mode shapes corresponding to each of the zero eigenvalues and are representative of a linear combination of the four pure mode shapes (x-translation, y-translation, rotational, and mechanism).

For the eigenvalue/eigenvector problem:

$$Kv = \lambda v$$

We can let $$\begin{Bmatrix} u_1 & u_2 & u_3 & u_4 \end{Bmatrix}$$ be the pure eigenvetors corresponding to the four zero eigenvalues. Assuming this tells us that:

$$Ku_i = 0u_i$$ for i = 1,2,3,4.

To get the true eigenvector, w, for each zero eigenvalue, we take the sum of the four contributors, each multiplied by a constant $$\alpha _{1}$$, giving us:

$$ w = \alpha_1 u_1 + \alpha_2 u_2 + \alpha_3 u_3 + \alpha_4 u_4$$

Because w consists of the combination of eigenvectors for zero eigenvalues, w too represents an eigenvector corresponding to a zero eigenvalue. That is that:

"$\bar{K}\bar{w} = K(\Sigma \alpha _iu_i) = \bar{0} = 0*\bar{w}$"

 Wednesday October 15 

Element Stiffness Matrix to Global Stiffness Matrix
Justification of assembly of element stiffness matrices into global stiffness matrix K.

$$K^{(e)} \; \; e=1, ... ,n_{el}\;$$

Consider the 2 bar truss example. Recall Force displacement equation.

$$ K^{(e)}_{4x4} \cdot d^{(e)}_{4x1} = F^{(e)}_{4x1}$$

You then can use the Euler cut method to establish the equilibrium of global node 2

Using the free body diagram of element one and of element 2 you can describe the elements degrees of freedom $$d^{(e)}_2 $$

$$d_2=d_2^{(1)}$$

$$d_3=d_3^{(1)}=d_1^{(2)} $$

$$d_4=d_4^{(1)}=d_2^{(3)}$$

$$d_5=d_3^{(2)}$$

$$d_6=d_4^{(2)}$$

Equilibrium of node 2





$$\sum{F_x}=-f_3^{(1)}-F_1^{(2)}=0$$

$$\sum{F_y}=P-f_4^{(1)}-F_2^{(2)}=0$$

Then you would use element force displacement relation $$K^{(e)}_{4x4} \cdot d^{(e)}_{4x1} = F^{(e)}_{4x1}$$

A Closer Look at the Force Displacement Relationship
Equation One:   $$\displaystyle F_3^{(1)}+F_1^{(2)}=0 $$

Equation Two:   $$ F_4^{(1)}+F_2^{(2)}=P $$

From the Force Displacement relationship: k(1)*d(1) = F(1), We can break the forces above into their k*d form:

$$ f_3^{(1)} = \begin{bmatrix} k_{31}^{(1)}d_1^{(1)} + k_{32}^{(1)}d_2^{(1)}+k_{33}^{(1)}d_3^{(1)}+k_{34}^{(1)}d_4^{(1)} \end{bmatrix} $$

$$ f_1^{(2)} = \begin{bmatrix} k_{11}^{(2)}d_1^{(2)} + k_{12}^{(2)}d_2^{(2)}+k_{13}^{(2)}d_3^{(2)}+k_{14}^{(2)}d_4^{(2)} \end{bmatrix} $$

Equation One Becomes:

$$ \begin{bmatrix} k_{31}^{(1)}d_1^{(1)} + k_{32}^{(1)}d_2^{(1)}+k_{33}^{(1)}d_3^{(1)}+k_{34}^{(1)}d_4^{(1)} \end{bmatrix} +  \begin{bmatrix} k_{11}^{(2)}d_1^{(2)} + k_{12}^{(2)}d_2^{(2)}+k_{13}^{(2)}d_3^{(2)}+k_{14}^{(2)}d_4^{(2)} \end{bmatrix} = 0 $$

And when we convert to Global dof's:

$$ \begin{bmatrix} k_{31}^{(1)}d_1 + k_{32}^{(1)}d_2+k_{33}^{(1)}d_3+k_{34}^{(1)}d_4 \end{bmatrix} +  \begin{bmatrix} k_{11}^{(2)}d_3 + k_{12}^{(2)}d_4 + k_{13}^{(2)}d_5 + k_{14}^{(2)}d_6 \end{bmatrix} = 0 $$

We can take the k coefficients of d3 and d4 and reduce the Matrices

This Proves our Global Stiffness K, from our earlier learning's of the Finite Element Method:

$$ K = \begin{bmatrix} k_{33}^{(1)}+k_{11}^{(2)} & k_{34}^{(1)}+k_{12}^{(2)} \\ k_{43}^{(1)}+k_{21}^{(2)} & k_{44}^{(1)}+k_{22}^{(2)} \end{bmatrix} $$

Another way to arrive at Global K
K= Ae * ke, where K is an n x n Matrix.

n is the total # of global dof's before eliminating the boundary conditions.

ned is the # of element dof's

Then:      Kn x n = A(e)undefined * k(e)ned x ned

A is the Assembly Operator

Virtual Work
The General Equation $$ P * \upsilon * \omega  $$

Define q:

q2x1 = T(e) * d(e)

Define k

k(e) = T(e) T * k(e) * T(e)

Deriving Finite Element Method PDE's
Equation (3) is Implied: Equation (3): $$ k * d - F = 0 $$

Equation (4) is derived and is not tirvial:

Equation (4): $$ w(k * d -F) = 0 $$, for all w

5 Truss MATLAB code
The 5 bar truss system in example 1.4 in the book can be solved by the below MATLAB code. analysis and explanation of the MATLAB code is shown below. This code and the the functions were downloaded from the books accompanying website:

http://bcs.wiley.com/he-bcs/Books?action=index&itemId=0471648086&bcsId=2256

Functions
The functions PlaneTrussElement.m, NodalSoln.m, PlaneTrussResults.m were all called in the code for the 5 truss system. They are all necessary for the completion of the truss problem.

PlaneTrussElement.m

This function take the Young's Modulus (e), the Area of the element (A) and the coordinates of the element ends (coord) and generates the stiffness matrix for a plane truss element. The function first calculates the lengths of the the elements. Once those are calculated the stiffness matrix is calculated.

NodalSoln.m

The function takes the global coefficient matrix (K), the global right hand side vector (R), the list of degrees of freedom with specified values, and the specified values to determine the displacements and reactions at each node. The dof of the the system is first found by using the  command which finds the longest dimension of R. df then finds the difference between the dofs that have known values (a value of zero) and the dof that were found in the previous line. The displacements and the reactions are then calculated.

PlaneTrussResults.m

This function computes the plane truss element results. It takes in the modulus of elasticity (e), the area of the cross-section (A), the coordinates at the element ends (coord), and the displacements at the elemtent ends (disps) and calculates the axial strain (eps), stress(sigma) and force (force). The results are stored in the variable  and sent back to the main program.

Results from the 5 Truss MATLAB example
The results from the 5 Truss System MATLAB code is shown below. As stated in the PlaneTrussResults.m function the  columns go in the order of axial strain, axial stress, and last axial force.

Code for the plot of the undeformed and deformed five truss system
Below is the code for the plot of the five bar system. Using the results that were obtained from the above MATLAB solution the deformed plot was created. The dotted lines are the undeformed system and the solid lines are the deformed system. The the displacements were multiplied by a factor of 500 so that they would appear on the plot.

Three Bar Truss System
A program was also written to solve the three bar truss system. The code follows the model code that was given to us in HW 2. The functions PlaneTrussElement.m, NodalSoln.m, and PlaneTrussResults.m were once again used in the code to obtain the results. Descriptions of these functions can be found in the section for the Five bar truss system.

Below are the results from the above code.

Plot of undeformed and deformed three bar truss system
Below is the code and the plot of the three bar truss system. The plot was made from the results that were found in the above code. The results demonstrate that only node two had a displacement of (-15.167, 48.231). The undeformed system is plotted with the dotted lines and the deformed system is plotted with solid lines. The deformed displacements have been scaled down by a tenth in order to enhance in the plot.

Contributing Team Members
Christopher Pergola EML4500.F08.RAMROD.B 18:10, 23 October 2008 (UTC)

Eml4500.f08.RAMROD.A 18:46, 23 October 2008 (UTC)

Chris Fontana Eml4500.f08.RAMROD.F 18:52, 24 October 2008 (UTC)

Eml4500.f08.ramrod.D 19:07, 23 October 2008 (UTC)

EML4500.f08.RAMROD.E 22:10, 23 October 2008 (UTC)

Eml4500.f08.ramrod.c 21:03, 24 October 2008 (UTC)