User:Eml5526.s11.team6.gravois

Team Homework Report Page Homework 1 Homework 2 Homework 3 Homework 4 = Homework 1 =

Find
Derive the following equation by performing a balance of forces:

$$ \frac{\partial}{\partial x}\left [ A(x)E(x)\frac{\partial u}{\partial x} \right ]+f(x,t)=A(x) \rho(x)\cdot \frac{\partial ^{2}u}{\partial t^{2}}\ $$

Solution
Consider the FBD of an infinitesimally small section of length dx of the 1-D elastic bar with varying material properties at a distance x from the origin. Where, A(x) - Varying cross section of the elastic bar, E(x) - Modulus of elasticity, rho (x) - Mass density, u - Displacement constrained only to x- axis (1-D bar), t - Time.


 * {| style="width:100%" border="0"

Where, $$\; {A}'=A(x) (dx)\cdot \rho (x)\cdot\frac{\partial^{2}u }{\partial t^{2}} =mass\cdot acceleration=Inertia\, Force$$
 * $$\displaystyle (Eq. 1.1) $$
 * }
 * }

As per Newton's law of motion, the involved forces should contribute to inertia force.

Force equilibrium in the x-direction yields:

$$\sum F_{x}=0 $$


 * {| style="width:100%" border="0"

$$\sum F_{x}=-N(x,t)+N(x+dx,t)+f(x,t) dx-m(x)\cdot \frac{\partial ^{2}u}{\partial t^{2}}=0$$
 * $$\displaystyle (Eq. 1.2) $$
 * }
 * }


 * {| style="width:100%" border="0"

Where, $$\;\; m(x)=A(x)(dx)\rho (x)$$
 * $$\displaystyle (Eq. 1.3) $$
 * }
 * }

Including Taylor series expansion for $$N(x+dx,t)$$ and neglecting the higher order terms (h.o.t) yields:

$$ \cancel {-N(x,t)}+[\cancel {N(x,t)}+\frac {\partial }{\partial x}(N(x,t))+(h.o.t)]+f(x,t)=\frac{m(x)}{dx}\cdot\frac{\partial ^{2} u}{\partial t^{2}}$$


 * {| style="width:100%" border="0"

$$\frac{\partial }{\partial x}(N(x,t))+f(x,t)=\frac{m(x)}{dx}\cdot\frac{\partial ^{2} u}{\partial t^{2}}$$
 * $$\displaystyle (Eq. 1.4) $$
 * }
 * }


 * {| style="width:100%" border="0"

$$Also,\; N(x,t)=A(x) \cdot \sigma (x)=A(x) \cdot [E(x)\cdot \epsilon (x,t)] $$
 * $$\displaystyle (Eq. 1.5) $$
 * }
 * }


 * {| style="width:100%" border="0"

$$\;\;\;\;\;\;\;\;\;\;\;\Rightarrow N(x,t)= A(x)\cdot E(x) \cdot\frac{\partial u}{\partial x}(x,t)$$
 * $$\displaystyle (Eq. 1.6) $$
 * }
 * }

Thus proving:


 * {| style="width:100%" border="0"




 * $$\displaystyle (Eq. 1.7) $$
 * }
 * }

Find
Discuss the case in which the bar has a rectangular cross section shown below.

Solution
We have from reference of class notes for mtg 2, eq(3) pg.2-3,


 * {| style="width:100%" border="0"

$$m=\rho (x) \cdot A(x) $$
 * $$\displaystyle (Eq. 2.1) $$
 * }
 * }

Also with reference to class notes for mtg 6, pg 6-1, point 2), for A(x)=bˑh(h),


 * {| style="width:100%" border="0"

$$m=\rho \left ( x+\frac{dx}{2} \right )\cdot \frac{1}{2}\left [ h(x)+h(x+dx) \right ]b $$
 * $$\displaystyle (Eq. 2.2) $$
 * }
 * }

Substituting Eq 2.1 and Eq 2.2 into Eq 1.7 for case of A(x)=bˑh(x),
 * {| style="width:100%" border="0"

$$\frac{\partial}{\partial x}\left [ \left ( b\cdot h(x)\cdot E(x)\cdot \frac{\partial u}{\partial x}\right ) \right ]+f(x,t)=\frac{1}{2}\left [ h(x)+h(x+dx) \right ]\rho (x+dx)\cdot \frac{\partial^{2} u}{\partial t^{2}}$$
 * $$\displaystyle (Eq. 2.3) $$
 * }
 * }

Note:

1) $$\left [ \frac{1}{2}\left [ h(x)+h(x+dx) \right ]\rho (x+dx) \right ]$$ denotes mass per unit length at a distance of (x+dx/2), which is the middle of the elementary section considered in Figure 1.2.

2) $$\frac {1}{2}\left [h(x)+h(x+dx) \right ]$$ denotes average height of elementary section of figure 1.2.

3) Since we are considering an infinitesimally small section (dx) of the entire bar, we can write,

$$h(x)\approx h(x+dx)$$

$$\Rightarrow \frac{1}{2}\left [ h(x)+h(x+dx) \right ]\approx h(x)$$

4) Also, b is the width of the cross section and is constant and ρ(x+dx) is approximately equal to ρ(x) for a small section as in 3) above.

Eq 2.3 converts to the following:


 * {| style="width:100%" border="0"




 * $$\displaystyle (Eq. 2.4) $$
 * }
 * }

Eq 2.4 is the modified form of Eq 1.7 for the specific case of A(x)=bˑh(x).

Eml5526.s11.team6.gravois 04:37, 25 January 2011 (UTC)

= Homework 2 =

Find
1) Construct Γ5x5(F) and observe the proportionality of Γ.

2) Find det Γ(F)

3) Conclude F is orthogonal basis i.e. $$\Gamma _{ij}=\left \langle b_{i},b_{j} \right \rangle=\delta _{ij}$$

Given
$$F=\left \{ 1,cos(\omega x),cos(2\omega x),sin(\omega x),sin(2\omega x) \right \} , \Omega =[0,T], i=1,2$$

Solution
1) $$ \mathbf{\Gamma} _{5x5}=\int_{\Omega }{}b_{i}\cdot b_{j}\;dx=\int_{0}^{T}b_{i}\cdot b_{j}\;dx $$

Constructing a matrix of bi ˑ bj yields:

$$ \Rightarrow \begin{bmatrix} 1\cdot 1 & cos(\omega x)\cdot 1 & cos(2\omega x)\cdot 1 & sin(\omega x)\cdot 1 & sin(2\omega x)\cdot 1\\ 1\cdot cos(\omega x) & cos(\omega x)\cdot cos(\omega x) & cos(2\omega x)\cdot cos(\omega x) & sin(\omega x)\cdot cos(\omega x) & sin(2\omega x)\cdot cos(\omega x)\\ 1\cdot cos(2\omega x) & cos(\omega x)\cdot cos(2\omega x) & cos(2\omega x)\cdot cos(2\omega x) & sin(\omega x)\cdot cos(2\omega x) & sin(2\omega x)\cdot cos(2\omega x)\\ 1\cdot sin(\omega x) & cos(\omega x)\cdot sin(\omega x) & cos(2\omega x)\cdot sin(\omega x) & sin(\omega x)\cdot sin(\omega x) & sin(2\omega x)\cdot sin(\omega x)\\ 1\cdot sin(2\omega x) & cos(\omega x)\cdot sin(2\omega x) & cos(2\omega x)\cdot sin(2\omega x) & sin(\omega x)\cdot sin(2\omega x) & sin(2\omega x)\cdot sin(2\omega x) \end{bmatrix} $$

Multiplying out the above matrix and integrating the result using a table of integrals yielded the following matrix:

$$ \Rightarrow \begin{bmatrix} x & sin(\omega x) & sin(\omega x)cos(\omega x) & -cos(\omega x) & -\frac{1}{2}cos(2\omega x)\\ sin(\omega x) & \frac{1}{2}(x+sin(\omega x)cos(\omega x)) & \frac{1}{6}(3sin(\omega x)+sin(3\omega x)) & -\frac{1}{2}cos^2(\omega x) & -\frac{2}{3}cos^3(\omega x)\\ sin(\omega x)cos(\omega x) & \frac{1}{6}(3sin(\omega x)+sin(3\omega x)) & \frac{1}{8}(4x+sin(4\omega x)) & \frac{1}{6}(3cos(\omega x)-cos(3\omega x)) & -\frac{1}{8}cos(4\omega x))\\ -cos(\omega x)) & -\frac{1}{2}cos^2(\omega x) & \frac{1}{6}(3cos(\omega x)-cos(3\omega x)) & \frac{1}{2}(x-sin(\omega x)cos(\omega x)) & \frac{2}{3}sin^3(\omega x)\\ -\frac{1}{2}cos(2\omega x) & -\frac{2}{3}cos^3(\omega x) & -\frac{1}{8}cos(4\omega x)) & \frac{2}{3}sin^3(\omega x) & \frac{1}{8}(4x-sin(4\omega x)) \end{bmatrix} $$

By definition of Angular Frequency

$$2\pi =\omega T\rightarrow \omega=\frac{2\pi}{T}$$

Solving the above matrix over the limits of the bar from 0 to T yielded:


 * {| style="width:100%" border="0"




 * }
 * }

Г is a symmetric matrix about the diagonal.

2) The determinate of the above Г matrix was computed using MATLAB as shown.




 * {| style="width:100%" border="0"




 * }
 * }

3) An orthogonal matrix is a matrix such that the transpose of the matrix multiplied by the original matrix should yield the identity matrix I, a matrix with only 1's on the diagonal and zeros everywhere else.




 * {| style="width:100%" border="0"




 * }
 * }

Find
1) Construct Γ5x5(F) and observe the proportionality of Γ.

2) Find det Γ(F)

3) Conclude F is orthogonal basis i.e. $$\Gamma _{ij}=\left \langle b_{i},b_{j} \right \rangle=\delta _{ij}$$

Given
$$F=\left \{ 1,x,x^2,x^3,x^4 \right \},

\Omega =[0,1]$$

Solution
1) $$ \mathbf{\Gamma} _{5x5}=\int_{\Omega }{}b_{i}\cdot b_{j}\;dx=\int_{0}^{1}b_{i}\cdot b_{j}\;dx $$

Constructing a matrix of bi ˑ bj yields:

$$ \Rightarrow \begin{bmatrix} 1\cdot 1 & x\cdot 1 & x^2\cdot 1 & x^3\cdot 1 & x^4\cdot 1\\ 1\cdot x & x\cdot x & x^2\cdot x & x^3\cdot x & x^4\cdot x\\ 1\cdot x^2 & x\cdot x^2 & x^2\cdot x^2 & x^3\cdot x^2 & x^4\cdot x^2\\ 1\cdot x^3 & x\cdot x^3 & x^2\cdot x^3 & x^3\cdot x^3 & x^4\cdot x^3\\ 1\cdot x^4 & x\cdot x^4 & x^2\cdot x^4 & x^3\cdot x^4 & x^4\cdot x^4 \end{bmatrix} $$

Multiplying out the above matrix yields:

$$ \Rightarrow \begin{bmatrix} 1 & x & x^2 & x^3 & x^4\\ x & x^2 & x^3 & x^4 & x^5\\ x^2 & x^3 & x^4 & x^5 & x^6\\ x^3 & x^4 & x^5 & x^6 & x^7\\ x^4 & x^5 & x^6 & x^7 & x^8 \end{bmatrix} $$

Integrating the above matrix with respect to x yields:

$$ \Rightarrow \begin{bmatrix} x & \frac{x^2}{2} & \frac{x^2}{2} & \frac{x^3}{3} & \frac{x^4}{4}\\ \frac{x^2}{2} & \frac{x^3}{3} & \frac{x^4}{4} & \frac{x^5}{5} & \frac{x^6}{6}\\ \frac{x^3}{3} & \frac{x^4}{4} & \frac{x^5}{5} & \frac{x^6}{6} & \frac{x^7}{7}\\ \frac{x^4}{4} & \frac{x^5}{5} & \frac{x^6}{6} & \frac{x^7}{7} & \frac{x^8}{8}\\ \frac{x^5}{5} & \frac{x^6}{6} & \frac{x^7}{7} & \frac{x^8}{8} & \frac{x^9}{9} \end{bmatrix}_{0}^{1} $$

Evaluating the integrated matrix with the limits of [0,1] yielded:

$$ \Rightarrow \begin{bmatrix} 1-0 & \frac{1^2}{2}-0 & \frac{1^3}{3}-0 & \frac{1^4}{4}-0 & \frac{1^5}{5}-0\\ \frac{1^2}{2}-0 & \frac{1^3}{3}-0 & \frac{1^4}{4}-0 & \frac{1^5}{5}-0& \frac{1^6}{6}-0\\ \frac{1^3}{3}-0 & \frac{1^4}{4}-0 & \frac{1^5}{5}-0 & \frac{1^6}{6}-0& \frac{1^7}{7}-0\\ \frac{1^4}{4}-0 & \frac{1^5}{5}-0 & \frac{1^6}{6}-0 & \frac{1^7}{7}-0& \frac{1^8}{8}-0\\ \frac{1^5}{5}-0 & \frac{1^6}{6}-0 & \frac{1^7}{7}-0 & \frac{1^8}{8}-0& \frac{1^9}{9}-0 \end{bmatrix} $$

The family of linearly independent basis functions matrix is thus:


 * {| style="width:100%" border="0"




 * }
 * }

Г is a symmetric matrix about the diagonal.

2) The determinate of the above Г matrix was computed using MATLAB as shown.




 * {| style="width:100%" border="0"




 * }
 * }

3) An orthogonal matrix is a matrix such that the transpose of the matrix multiplied by the original matrix should yield the identity matrix I, a matrix with only 1's on the diagonal and zeros everywhere else.




 * {| style="width:100%" border="0"




 * }
 * }

Eml5526.s11.team6.gravois 01:01, 1 February 2011 (UTC)

=Homework 3=

Given

 * {| style="width:100%" border="0"


 * $$\alpha =b_i \left [ a_{2}'b_{j}'+a_{2}b_{j}''\right ]$$
 * $$\displaystyle (Eq. 4.1) $$
 * }
 * {| style="width:100%" border="0"
 * }
 * {| style="width:100%" border="0"

Where,
 * $$\beta =b_j \left [ a_{2}'b_{i}'+a_{2}b_{i}''\right ]$$
 * $$\displaystyle (Eq. 4.2) $$
 * }
 * }
 * }


 * {| style="width:100%" border="0"

$$\overline{K}_{ij}=\int_{\Omega }b_i (x)\left [ \frac{d}{dx}\left (a_2\frac{d}{dx}b_j (x) \right ) \right ]dx$$
 * $$\displaystyle (Eq. 4.3) $$
 * }
 * }


 * {| style="width:100%" border="0"

$$\overline{K}_{ji}=\int_{\Omega }b_j (x)\left [ \frac{d}{dx}\left (a_2\frac{d}{dx}b_i (x) \right ) \right ]dx$$
 * $$\displaystyle (Eq. 4.4) $$
 * }
 * }

α represents the integrand in the stiffness matrix (Eq 4.3) and β represents the integrand in stiffness matrix (Eq 4.4).

Find
Show that the above $$\alpha \neq \beta $$ using the following example:

$$\left\{\begin{matrix} b_i(x)=cos(ix)\\ b_j(x)=cos(jx) \end{matrix}\right.$$ where $$i \neq j$$

Solution
Substituting bi and bj into the α and β equations yielded:

$$\alpha = cos(ix)\left [ a_2'sin(jx)'+a_2cos(jx)'' \right ]$$ $$\beta = cos(jx)\left [ a_2'sin(ix)'+a_2cos(ix)'' \right ]$$

Carrying out the derivations yielded the following:

$$\alpha = cos(ix)\left [ -ja_2'sin(jx)-j^{2}a_2cos(jx) \right ]$$ $$\beta = cos(jx)\left [ -ia_2'sin(ix)-j^{2}a_2cos(ix) \right ]$$

Thus, $$\alpha \neq \beta$$

Since $$\alpha \neq \beta$$, it proves that the stiffness matrix $$\overline{K}$$ is not symmetric.

Given
The three-bar structure below is subjected to the given load of 103N at point B. Th Young's modulus is E=1011Pa, the cross sectional area of the bar BC is 2x10-2m2 and BD is 10-2m2. Point D is free to move in the x-direction and all distances are given in meters.



Find
a. Construct the global stiffness matrix and load matrix. b. Partition the matrices and solve for the unknown displacements at Point B and Point D. c. Find the stresses in the three bars. d. Find the reactions at nodes C, D, F.

a. Construct the global stiffness matrix and load matrix.
The stiffness matrix for each bar element is given below:

$$\mathbf{K}^{(n)}=\frac{AE}{l}\begin{bmatrix} cos(\phi)^2 & cos(\phi)sin(\phi) & -cos(\phi)^2 & -cos(\phi)sin(\phi)\\ cos(\phi)sin(\phi) & sin(\phi)^2 & -cos(\phi)sin(\phi) & -sin(\phi)^2\\ -cos(\phi)^2 & -cos(\phi)sin(\phi) & cos(\phi)^2 & cos(\phi)sin(\phi)\\ -cos(\phi)sin(\phi) & -sin(\phi)^2 & cos(\phi)sin(\phi) & sin(\phi)^2 \end{bmatrix}$$



For bar (1):

$$cos(\phi)=cos(90)=0$$ $$sin(\phi)=sin(90)=1$$

Substituting the above results into the general element stiffness matrix yielded the following:

$$\mathbf{K}^{(1)}=\frac{2AE}{l}\begin{bmatrix} 0 & 0 & 0 & 0\\ 0 & 1 & 0 & -1\\ 0 & 0 & 0 & 0\\ 0 & -1 & 0 & 1 \end{bmatrix}$$

For bar (2):

$$cos(\phi)=cos(45)=\frac{1}{\sqrt{2}}$$ $$sin(\phi)=sin(45)=\frac{1}{\sqrt{2}}$$

Substituting the above results into the general element stiffness matrix yielded the following:

$$\mathbf{K}^{(2)}=\frac{AE}{\sqrt{2}l}\begin{bmatrix} \frac{1}{2} & \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2}\\ \frac{1}{2} & \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2}\\ -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2}\\ -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \end{bmatrix}$$

For bar (3):

$$cos(\phi)=cos(135)=-\frac{1}{\sqrt{2}}$$ $$sin(\phi)=sin(135)=\frac{1}{\sqrt{2}}$$

Substituting the above results into the general element stiffness matrix yielded the following:

$$\mathbf{K}^{(3)}=\frac{AE}{\sqrt{2}l}\begin{bmatrix} \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2}\\ -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} & -\frac{1}{2}\\ -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} & -\frac{1}{2}\\ \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2} \end{bmatrix}$$

Inserting the above element stiffness matrices into the global stiffness matrix yielded the following:

$$\mathbf{K}=\frac{AE}{l}\begin{bmatrix} 0+\frac{1}{2\sqrt{2}}+\frac{1}{2\sqrt{2}} & 0+\frac{1}{2\sqrt{2}}+\frac{1}{2\sqrt{2}} & 0 & 0 & -\frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}}\\ 0+\frac{1}{2\sqrt{2}}-\frac{1}{2\sqrt{2}} & 2+\frac{1}{2\sqrt{2}}+\frac{1}{2\sqrt{2}} & 0 & -2 & \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}}\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & -2 & 0 & 2 & 0 & 0 & 0 & 0\\ \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & 0 & 0\\ \frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & 0 & 0\\ -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & \frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}}\\ -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & \frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} \end{bmatrix}$$

After combining terms in the global stiffness matrix, the result is as follows:


 * {| style="width:100%" border="0"




 * }
 * }

The general load matrix is given below:

$$\mathbf{F}=\begin{bmatrix} F_{1x} & F_{1y} & F_{2x} & F_{2y} & F_{3x} & F_{3y} & F_{4x} & F_{4y}\\ \end{bmatrix}^{T}$$

Substituting in the force acting at Point B yielded the load matrix below:


 * {| style="width:100%" border="0"




 * }
 * }

b. Partition the matrices and solve for the unknown displacements at Point B and Point D.
The Global systems of equations become:

F = K d

$$\begin{bmatrix} 1000 \\ 0 \\ r_{2x} \\ r_{2y} \\ r_{3x} \\ r_{3y} \\ 0 \\ r_4y\\ \end{bmatrix}=1x10^9\begin{bmatrix} \frac{1}{\sqrt{2}} & 0 & 0 & 0 & -\frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}}\\ 0 & 2+\frac{1}{\sqrt{2}} & 0 & -2 & \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}}\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & -2 & 0 & 2 & 0 & 0 & 0 & 0\\ \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & 0 & 0\\ \frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & 0 & 0\\ -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & \frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}}\\ -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & 0 & 0 & 0 & 0 & \frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} \end{bmatrix}\cdot\begin{bmatrix} u_{1x}\\ u_{1y}\\ 0\\ 0\\ 0\\ 0\\ u_{4x}\\ 0

\end{bmatrix}$$

Partitioning to solve for the displacement of Point B yields the following matrices:

$$\begin{bmatrix} 1000\\0 \\0 \end{bmatrix}=10^9 \begin{bmatrix} \frac{1}{\sqrt2} & 0 & -\frac{1}{2\sqrt2}\\ 0 & 2+\frac{1}{\sqrt2} & -\frac{1}{2\sqrt2}\\ -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \end{bmatrix}\cdot \begin{bmatrix} u_{1x}\\u_{1y} \\u_{4x} \end{bmatrix}$$

$$\begin{bmatrix} 1000\\0 \\0 \end{bmatrix}= \begin{bmatrix} 7.071x10^8 & 0 & -3.536x10^8\\ 0 & 2.7071x10^9 & -3.536x10^8\\ -3.536x10^8 & -3.536x10^8 & 3.536x10^8 \end{bmatrix}\cdot \begin{bmatrix} u_{1x}\\u_{1y} \\u_{4x} \end{bmatrix}$$

Solving for the displacements yields:

$$\mathbf{d}=\begin{bmatrix} 3.328x10^{-6} \\ 5x10^{-7} \\ 3.828x10^{-6} \end{bmatrix}m$$


 * {| style="width:100%" border="0"




 * }
 * }

c. Find the stresses in the three bars.
The stress in each member is obtained using the following equation:

$$\sigma^{(1)}=\frac{E^{(n)}}{l^{(n)}}\begin{bmatrix} -cos(\phi^{(n)})& -sin(\phi^{(n)})& cos(\phi^{(n)})& sin(\phi^{(n)}) \end{bmatrix}\mathbf{d}^{(n)}$$

For Bar 1:

$$\sigma^{(1)}=\frac{10^9}{1}\begin{bmatrix} 0 & -1& 0& 1 \end{bmatrix}\mathbf\cdot \begin{bmatrix} 3.328x10^{-6}\\ 5x10^{-7}\\ 0\\ 0 \end{bmatrix}$$

For Bar 2:

$$\sigma^{(2)}=\frac{10^9}{1}\begin{bmatrix} -\frac{1}{\sqrt2} & -\frac{1}{\sqrt2}& \frac{1}{\sqrt2}& \frac{1}{\sqrt2} \end{bmatrix}\mathbf\cdot \begin{bmatrix} 3.328x10^{-6}\\ 5x10^{-7}\\ 3.828x10^{-6}\\ 0 \end{bmatrix}$$

For Bar 3:

$$\sigma^{(2)}=\frac{10^9}{1}\begin{bmatrix} \frac{1}{\sqrt2} & -\frac{1}{\sqrt2}& -\frac{1}{\sqrt2}& \frac{1}{\sqrt2} \end{bmatrix}\mathbf\cdot \begin{bmatrix} 3.328x10^{-6}\\ 5x10^{-7}\\ 0\\ 0 \end{bmatrix}$$

d. Find the reactions at nodes C, D, F.
$$\mathbf{r}=\begin{bmatrix} r_{2x}\\ r_{2y}\\ r_{3x}\\ r_{3y}\\ r_{4y}\\ \end{bmatrix}=\mathbf{Kd}=10^{9}\begin{bmatrix} 0 &0 & 0\\ 0 & -2 & 0\\ -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2}\\ -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0\\ \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 \end{bmatrix}\cdot \begin{bmatrix} 3.328x10^{-6}\\ 0.5x10^{-6}\\ 3.828x10^{-6} \end{bmatrix}= \begin{bmatrix} 0\\ -1000\\ 0\\ -1000\\ 1000 \end{bmatrix}N $$

= Homework 4 =

Given
Consider the elastic bar with a variable distributed spring p(x) along its length shown below in Figure 4.6.1. The distributed distributed s[pring imposes an axial force on the bar in proportion to the displacement. The bar is of length l, cross-sectional area A(x), Young's modulus E(x) with body force b(x) and boundary conditions shown in the following figure.



The above problem statement is given in A First Course in Finite Elements FB p74 Problem 3.9.

Find
a. Construct the strong form. b. Construct the weak form.

Solution
The bar must satisfies the following conditions: 1. The bar must be in equilibrium. 2. The bar must satisfy Hooke's Law σ(x)=E(x)ε(x). 3. The displacement field must be compatible. 4. It must satisfy the strain-displacement equation.

The following Figure 4.6.2 shows the internal force balance of a segment of the bar.



a. Strong Form
Performing a force balance on the section in Figure 4.6.2 yields the following differential equation 4.6.1.

Dividing the above equation through by dx to simplify terms yielded the following equation 4.6.2.

Now take the limit as Δx goes to 0. The result is presented below in equation 4.6.3.

By Hooke's Law:

Substituting Eq.4.6.6 into Eq.4.6.5 and then plugging that result into Eq.4.6.4 yielded the following equation.

Substituting Eq.4.6.7 into Eq.4.6.3 yields the final strong for differential equation as follows:


 * {| style="width:100%" border="0"




 * }
 * }

b. Weak Form
For the weak form, we multiply the strong form solution by the weighting function w and integrate the result over the length of the bar 0 to l as shown in the following equation:

Integrating the first term of Eq.4.6.8 by using integration by parts yielded the following equation:

Knowing that: $$A(x)E(x)\frac{du}{dx}=\bar{t}$$ and accounting for the b(x)-p(x) term yielded the following:

The final weak form becomes the following:


 * {| style="width:100%" border="0"




 * }
 * }

Figure
Following figure 4.7.8 shows the geometry and the mesh generated by CGX for a sphere volume segment



Geometry generation code
Following figure 4.7.10 shows the code to generate the geometry and the mesh by CGX for the cylinder in Figure 4.7.9.

Explanation of the commands used is given in the succeeding subsection.



Explanation of the commands used
Please refer to the 'Disc' section above for explanation of all the following CGX commands used to build the geometry and mesh of the cylinder Following commands are used to generate the geometry and mesh of the cylinder in figure 4.7.3 -

' PNT  For entering a point - (refer to the 'Disc' section for explanation)

' LINE  For entering a Line - (refer to the 'Disc' section for explanation)

' GSUR  For entering a surface - (refer to the 'Disc' section for explanation)

' ELTY  To assign a specific element type to a set of entities and generating the mesh - (refer to the 'Disc' section for explanation)

=Homework 5=

Given
The strong form of GM1.0b/d1b on page 2 of Mtg 9 shown below

$$ \frac{\partial}{\partial x} \left[ (2+3x) \frac{\partial u}{\partial x} \right] + 5x = 0 $$,$$ \forall x \in \left[ 0,1 \right]$$, and $$\Gamma_g = [0]$$, is $$u(0)=4$$, $$\Gamma_h = [1]$$, is $$ -\frac{\partial u(x=1)}{\partial x}=6 $$

Solution
Using weighting function $$ w(x) $$ such as $$ w(0) = 0 $$ and the technique of integration by parts, we have derived the following weak form.

Fourier basis function
The trial solution then becomes the following for n=1 to convergence:

The essential boundary condition,$$u(x=0)=4$$ must be satisfied, and the basis functions must be selected so $$ w(n=0) = 0 $$. Solving the above basis function for u(x=1)=4 yielded the following: $$ {d}_{1}=4, i=1,2,\cdots $$, $$ j=1,2,\cdots $$

The remaining displacement coefficients must be solved for $$d_{j}, \;\;j=1,2,\cdots $$, by constructing the following stiffness and force matrices using the switching Fourier basis function.

Matlab was used to generate the following stiffness and force matrices shown below.

Since Kd=F, d can be solved for using Matlab to perform the matrix algebra. The resulting d matrix is shown below.

The trial solution of the governing equation then becomes the following.

Now we plot for Exact solution u(x) and Approximate solution uh(x) and also Convergence of error vs. number of degrees of freedom (n)

As seen from the Figure 5.1.1 below, it is evident that the approximate solution uh(x) very closely imitates the exact solution u(x) Also it can be seen from Figure 5.1.2. that the error converges to almost zero value when more than 8 basis functions are used





Find
Plot three figures similar to those of the linear Lagrangian element basis functions..

Solution
The generic equation for the quadratic Lagrangian element basis functions (QLEBF) is presented below in Equation 5.6.1;


 * {| style="width:100%" border="0"

$$ L_{j,n}=\prod_{k=1}^{n}\frac{(x-x_{k})}{(x_{i}-x_{k})} $$
 * $$\displaystyle (Eq.5.6.1) $$
 * }
 * }

Equation 5.6.1 can be expanded to the following equation to fit each of our cases as shown in Equations 5.6.2-4 below with x1=1, x2=2, and x3=3.


 * {| style="width:100%" border="0"

$$ L_{1,3}=\frac{(x-x_{2})(x-x_{3})}{(x_{1}-x_{2})(x_{1}-x_{3})} $$
 * <p style="text-align:right;">$$\displaystyle (Eq.5.6.2) $$
 * }
 * }


 * {| style="width:100%" border="0"

$$ L_{2,3}=\frac{(x-x_{3})(x-x_{1})}{(x_{2}-x_{3})(x_{2}-x_{1})} $$
 * <p style="text-align:right;">$$\displaystyle (Eq.5.6.3) $$
 * }
 * }


 * {| style="width:100%" border="0"

$$ L_{3,3}=\frac{(x-x_{1})(x-x_{2})}{(x_{3}-x_{1})(x_{3}-x_{2})} $$
 * <p style="text-align:right;">$$\displaystyle (Eq.5.6.4) $$
 * }
 * }

Since the quadratic Lagrangian element basis functions have three nodes, three QLEBFs were used to generate the following plots 5.6.1-4 in MATLAB. The code used to generate these plots is presented below.









Comments
Each of the plotted quadratic Lagrangian element basis functions has three nodes as opposed to the two nodes of the linear basis functions due to the n=3 condition. It is observed from the plots that the sum of the basis functions at each node must be equal to 1 regardless of the number of nodes. This also holds true for the linear Lagrangian element basis functions as well.

=Homework 6=

Given
The G2Dm1.0/D1 PDE

The G2Dm1.0/D1 PDE can be reduced to where $$\left[ \right]$$ is an identity matrix.

The essential boundary condition is $$g = 2\ \ \ \ on\ \ \ \ \partial \Omega $$.

Find
1) Solve G2DM1.0/D1 using 2D LIBF until accuracy $${10^{ - 6}}$$ at center (x,y)=(0,0) and verify results with any FE codes with detailed documentation.

2) Plot resulting solution in 3d with contours.

Solution
1) For Eq.6.5.1, we can write the weak form

And then we choose Lagrange Interpolation Basis Function for each element with n=m=2


 * $$N_{I}(x,y)=L_{i,m}(x)\cdot L_{j,n}(y),\;\; where \;\;I=i+(j-1)m$$

After substituting Eq.6.5.4 into Eq.6.5.3 in which Lagrange Interpolation Basis Function are both trial function and weight function, we integrate over the whole domain to get stiffness matrix and achieve final results using the following Matlab codes

The resulting stiffness and force matrices for the n=m=2 case is as follows:

$$K_{EF}=\begin{bmatrix} 0.6667 & -0.1667 & 0 & -0.1667 & 0 & 0 & 0 & 0 & -0.3333\\ -0.1667 & 1.3333 & -0.1667 & -0.3333 & -0.3333 & 0 & 0 & 0 & -0.3333\\ 0 & -0.1667 & 0.6667 & 0 & -0.1667 & 0 & 0 & 0 & -0.3333\\ -0.1667 & -0.3333 & 0 & 1.3333 & 0 & -0.1667 & -0.3333 & 0 & -0.3333\\ 0 & -0.3333 & -0.1667 & 0 & 1.3333 & 0 & -0.3333 & -0.16667 & -0.3333\\ 0 & 0 & 0 & -0.1667 & 0 & 0.6667 & -0.16667 & 0 & -0.3333\\ 0& 0 & 0 & -0.3333 & -0.3333 & -0.1667 & 1.3333 & -0.16667 & -0.3333\\ 0 & 0 & 0 & 0 & -0.1667 & 0 & -0.1667 & 0.6667 & -0.3333\\ -0.3333 & -0.3333 & -0.3333 & -0.3333 & -0.3333 & -0.3333 & -0.3333 & -0.3333 & 2.6667 \end{bmatrix}$$

$$F=\begin{bmatrix} 0.1110x10^{-15}& -0.2220x10^{-15}& 0.1110x10^{-15}& 0x10^{-15}& -0.1110x10^{-15}& 0.1110x10^{-15}& 0.1110x10^{-15}& 0x10^{-15}& \end{bmatrix}^{T}$$

Knowing that Kd=F and solving for the d matrix yielded:


 * $$d=

\begin{bmatrix} 2& 2& 2& 2& 2& 2& 2& 2& 2& \end{bmatrix}^{T}$$

The result from the MATLAB code is as follows:

The following is the matlab code used to generate the following plots as well as calculate the above matrices.

2)



Based on the above results we see that regardless of how many elements we use, we get a uniform temperature distribution. The exact solution is easily obtained, because this is static problem with uniform boundary conditions. The 3D contour plots are flat surfaces with no contour due to the uniform temperature distribution all the way around the perimeter of the plate. If the temperature varied along one of the edges or at a node, then the plots would have some contour and an uneven temperature distribution.

=My Team Members= Nakul Deshpande Shreyas Joglekar William Kurth Yongtao Li Anand Tupsakhare Chris Vork