text:
User:EAS4200C.F08.WIKI.A/Homework6

=Equation of Equilibrium Continued=

Recall the 1-D case from earlier: dimensional analysis.

$$\left. [f]=\frac{F}{L} \right \} \;\left [ \frac{f}{A} \right] = \frac{f}{L^3}$$

Notes :

 * The square brackets mean "dimension of"

 * F = Force

 * L = Length

$$[A]=L^2$$

$$[\sigma] = \frac{F}{L^2} \Rightarrow \left [ \frac{d \sigma}{dx} \right] = \frac{[d \sigma]}{[dx]} = \frac{\frac{F}{L^2}}{L} = \frac{F}{L^3}$$
 * F/L3 = Force/Volume or Body Force

$$[dx] = L$$

$$\left. \begin{matrix} \varepsilon = \frac{du}{dx}\\ E = \frac{\Delta L}{L} \end{matrix}  \right \} [\varepsilon] = \frac{[du]}{[dx]} = \frac{L}{L} = 1$$

$$\nu =\frac{\varepsilon_{yy}}{\varepsilon_{xx}} \Rightarrow [\nu] = \frac{[\varepsilon_{yy}]}{[\varepsilon_{xx}]} = 1$$

=Torsional Analysis Continued=

The torsional analysis previously started is continued by looking at the dimensions of the derivative of stress:

$$[\frac{\partial \sigma _{ij}}{\partial x_{i}}] = \frac{F}{L^3} = \frac{force}{volume}$$

Recalling - D. Prandtl stress function - from the road map previously constructed in class we see that:

$$\sigma _{yx} = \frac{\partial \phi }{\partial z}$$

and

 $$ \sigma _{zx} = -\frac{\partial \phi }{\partial y}$$

where $$\phi$$ plays the role of a potential function and $$(\sigma _{yx}, \sigma _{zx}$$ are composed of the "gradient" of $$\phi$$ with respect to (y,z).

As a review of the gradient operators, recall that for scalar functions f(x,y,z):  $$\bigtriangledown f(x,y,z) = \frac{\partial f}{\partial x}\vec{i} +\frac{\partial f}{\partial y}\vec{j} +\frac{\partial f}{\partial z}\vec{k}$$

Where f is the potential function.

Now, substituting the values for $$\sigma _{yx}$$ and $$\sigma _{zx}$$ into the previously given Laplace's Equation we get:

$$\frac{\partial }{\partial y}(\frac{\partial \phi}{\partial z}) + \frac{\partial }{\partial z}(-\frac{\partial \phi}{\partial y}) = 0$$

Which can be rearranged to show that:  $$\frac{\partial ^2 \phi}{\partial y \partial x} = \frac{\partial ^2 \phi}{\partial z \partial y}$$

Due to $$\phi$$ being a continuous and smooth function, the second derivatives are interchangeable in the order they are performed.

Likewise:

$$\frac{\partial ^2 \phi}{\partial y^2} + \frac{\partial ^2 \phi}{\partial z^2} = -2G\theta $$

The lefthand side can also be written as $$\bigtriangledown ^2 \phi$$ which is the Laplacian of $$\phi$$.

$$[t]_{3x1} = [\sigma]_{3x3}[n]_{3x1}$$

Where t is the component of the traction force t, $$\sigma$$ is the component of the stress tensor, and n is the component of the normal vector $$\vec{n}$$

A graphical depiction of the 2-D case is shown below:

2-D Case Continued

Here we will find the sum of the forces in the y-direction.
$$\displaystyle \begin{Vmatrix} \bar{n} \end{Vmatrix}= 1$$ $$\displaystyle \sum{F_y} = -\sigma_{yy}\cdot (dz\cdot 1) - \sigma(dy\cdot1) + t_y(ds\cdot1) = 0 = \begin{cases}dz =ds\cos \theta \\ dy = ds\sin\theta \end{cases}$$

$$\displaystyle 0 = -\sigma_{yy}n_y - \sigma_yzn_z + t_y ds$$

$$\displaystyle n_y = \cos\theta, n_z = \sin\theta $$ is the unit normal vector unit depth along x-direction.

This is equation one which will be combined with two later equations to give us the generalized 3-D case:

$$t_y=\sigma_{yy}n_y+\sigma_{yz}n_z\;$$

Note: $$\left[t_y \right]=\frac{F}{L^2}$$

$$\vec{t}\;$$:traction vector (distribution surface force)

$$\left[t_y \right]=\left[\sigma \right]$$

This is equation 2:

Taking the sum of the forces in the z-direction

$$t_z=\sigma _{yz}\cdot n_y+\sigma _{zz}n_z$$

This is equations one and two combined into matrix format:
 $$\begin{Bmatrix} t_y\\ t_z \end{Bmatrix}=\begin{bmatrix} \sigma_{yy} &\sigma_{yz} \\ \sigma_{zy}& \sigma_{zz} \end{bmatrix}\begin{Bmatrix} n_y\\ n_z \end{Bmatrix}$$

Generalized to 3-D case

$$\begin{Bmatrix} t_1\\ t_2\\ t_3 \end{Bmatrix}=\begin{bmatrix} \sigma_{11} & \sigma_{12} & \sigma_{13}\\ \sigma_{21} &\sigma_{22} & \sigma_{23}\\ \sigma_{31} & \sigma_{32} & \sigma_{33} \end{bmatrix} \begin{Bmatrix} n_1\\ n_2\\ n_3 \end{Bmatrix}$$

 $$t_i=\sum_{j=1}^{3}{\sigma_{ij}n_j}$$ $$i=1,2,3$$

 $$\left\{ti \right\}_{3x1}=\left[\sigma_{ij} \right]_{3x3}\left\{n_j \right\}_{3x1}$$

Two Cases of Torsional Analysis
There are two specific cases that we will concern ourselves with. The first being the closed thin walled cross-section and the second being the solid cross section. the two are depicted below.

Case 1: Closed Thin Walled Cross-Section (compare equations)
This case can be compared to the NACA Airfoil that is continuously being analyzed. It is important to note that the Prandtl stress function $$\displaystyle \Phi$$ is constant on the outer surface $$\displaystyle S_o$$ and the inner surface $$\displaystyle S_i$$. $$\displaystyle S_i$$ can also be represented by $$\displaystyle S_1$$ as well. The image below represents this case. It is comparable to Figure 3.11 in the text (Sun 2006).

Case 2: Solid Cross-Section
The one thing that is important to note about this case is that the Prandtl stress function $$\displaystyle \Phi$$ is zero on the outter surface $$\displaystyle S_o$$. This case is depicted below.

Uniform Bar With Solid Circular Cross-Section
This section coincides with section 3.3 in the text (Sun 2006). The following figure depicts the cross section of a solid circular bar and we will be determining the expression for the torque and polar moment of inertia in terms of the the Prandtl stress function.

Torque In Terms of $$\displaystyle \phi$$
The axis of the bar coincides with the origin of the coordinate system. The boundary of the bar can be expressed as the following equation;

$$\displaystyle x^2 +y^2 = a^2$$ where $$\displaystyle a^2$$ is the radius of the cross section

It is important to note that this is equivalent to the equation of a circle. We will assume that the stress function equals the following;

 $$\displaystyle \phi = C \left(\frac{x^2}{a^2} + \frac{y^2}{a^2} -1 \right)$$

With this information we are able to compute the torque. This is done below.

$$\displaystyle T = 2\int \int_{A}^{}{\phi dxdy} =2 C \int \int_{A}^{}{\left(\frac{x^2}{a^2} + \frac{y^2}{a^2} -1 \right) dxdy} = 2C\int \int_{A}^{}{\left(\frac{r^2}{a^2} -1 \right) dA} = 2C \left(\frac{{\pi} a^2}{2} - ({\pi}a^2) \right) = -C \pi a^2$$

We must now solve for the variable $$\displaystyle C$$. We know that the Prandtl stress function is equal to zero on the outer surface. We also that the gradient of the Prandtl stress function is equal to the following;

$$\displaystyle \nabla^2 \phi = -2G\theta = \frac{\partial^2 \phi }{\partial x^2} + \frac{\partial^2 \phi}{\partial y^2} = C\frac{2}{a^2} + C\frac{2}{a^2} = C\frac{4}{a^2}$$

Solving for $$\displaystyle C$$ yields;

 $$\displaystyle C = -\frac{a^2G\theta}{2}$$

Plugging this back into the equation for torque we can obtain the final expression.

$$\displaystyle T = \frac {\theta G \pi a^4}{2}$$

Polar Moment of Inertia In Terms of $$\displaystyle \phi$$
We know that the polar moment of inertia $$\displaystyle J$$ equals the following;

 $$\displaystyle J = \frac {T}{G\theta}$$

where

$$\displaystyle T = 2\int \int_{A}^{}{\phi dxdy}= -C\pi a^2$$

$$\displaystyle G \theta = -\frac{\nabla^2 \phi}{2} = \frac{2C}{a^2}$$

The polar moment of inertia then equals the following;

$$\displaystyle J = \frac {T}{G\theta} = --\frac{C\pi a^2}{\frac{2C}{a^2}} = \frac{\pi a^4}{2}$$

This is the exact equation that has been derived multiple time before using other methods.

{| class="toccolours collapsible collapsed" style="width:100%" ! Relationship Between K and k_c
 * The function K can be expressed in terms of k_c as the following;
 * The function K can be expressed in terms of k_c as the following;

$$\displaystyle K = k_c \frac{\pi^2}{12}$$

How do we come to this equation? First we have to recognize that the critical buckling load is given as follows;

$$\displaystyle (P_{x})_{cr} = k_c\frac{\pi^2D}{b}$$ where b is the width of the plate

The variable $$\displaystyle k_c = \left(\frac {mb}{a}+\frac{a}{mb} \right)^2$$. $$\displaystyle m $$ and $$\displaystyle n $$ are the number of half wave lengths and take on the value of whole numbers only. We next have to note that the critical buckling stress is given by the follwing equation;

$$\displaystyle (\sigma_{xx})_{cr} = \frac{(P_x)_{cr}}{bh} = k_c \frac{\pi^2D}{b^2h} $$ D is known as the plate bending stiffness

$$\displaystyle D = \frac{Eh^3}{12(1-v^2)}$$

For the case of a simply supported rectangular plate under a distributed compressive load on both sides, the critical buckling stress can be represented as the following;

$$\displaystyle (\sigma_{xx})_{cr} = K\frac{E}{1-v^2}\left(\frac{t}{b} \right)^2$$ where E is the elastic modulus, v is poissons ratio, and t is the thickness of the plate h.

Equating the two expressions for the critical buckling stress and solving for K, and making sure to plug in for the plate bending stiffness D, the following results can be obtained. As shown most of the variables cancel and

$$\displaystyle K\frac{E}{1-v^2}\left(\frac{h}{b} \right)^2 = k_c \frac{\pi^2Eh^3}{12b^2h(1-v^2)}$$

 $$\displaystyle K = k_c\frac{\pi^2}{12}$$

Verification
We will verify $$\displaystyle K $$ for the following values:
 $$\displaystyle \frac{a}{b} = 1, m=1$$

The value of $$\displaystyle k_c $$ would then equal the following;

$$\displaystyle k_c = \left(\frac {mb}{a}+\frac{a}{mb} \right)^2 = \left(\frac {(1)(1)}{1}+\frac{1}{(1)(1)} \right)^2 = 4$$

Plugging this value into the equation for $$\displaystyle K $$ we can obtain the final answer.

$$\displaystyle K = k_c\frac{\pi^2}{12} = 4\frac{\pi^2}{12} = 1.0472$$
 * }

Torque can be expressed as

 $$\displaystyle T = 2 \int_A \phi dA $$

From the theory of elasticity, the torque becomes

 $$\displaystyle T = 2C ( {J \over a^2} - A) $$

Where the area, $$\displaystyle A $$, since the radius is $$\displaystyle a $$, is simply

 $$\displaystyle A = \pi a^2 $$

And where $$\displaystyle J $$ is

 $$\displaystyle J = \int_A r^2 dA = {1 \over 2} \pi a^4 $$

Using the ad-hoc method

 $$\displaystyle T = GJ \theta $$

Recall that the Prandtl Stress function of $$\displaystyle \phi $$ states that the stress in the x-direction normal to the y-axis is

 $$\displaystyle \sigma_{yx} = {\partial \phi \over \partial z} $$

Taking this partial derivative, you obtain

 $$\displaystyle \sigma_{yx} = 2C {z \over a^2} $$

Where substitution for $$\displaystyle C $$ yields

 $$\displaystyle \sigma_{yx} = - G \theta z $$  (1.1)

Similarly, the stress in the x-direction normal to the z-axis is

 $$\displaystyle \sigma_{zx} = - {\partial \phi \over \partial y} $$

Again, taking this partial derivative, you obtain

 $$\displaystyle \sigma_{zx} = -2C {y \over a^2} $$

Where, again, substitution for $$\displaystyle C $$ yields

 $$\displaystyle \sigma_{zx} = G \theta y $$  (1.2)

Also note that shear stress, $$\displaystyle \tau $$, is

 $$\displaystyle \tau = {Tr \over J} $$

Proof of No Warping
Now, using the above equations in addition to the stress strain relation, we can prove that no warping takes place.

First, recall that

 $$\displaystyle \gamma = {\sigma \over G} $$

Or more specifically in the context of this problem

 $$\displaystyle \gamma_{yx} = {\sigma_{yx} \over G} $$

 $$\displaystyle \gamma_{zx} = {\sigma_{zx} \over G} $$

Also recall that

 $$\displaystyle \gamma_{yx} = {\partial u_x \over \partial y} - \theta z $$

 $$\displaystyle \gamma_{zx} = {\partial u_x \over \partial z} + \theta y $$

Thus,

 $$\displaystyle {\sigma_{yx} \over G} = {\partial u_x \over \partial y} - \theta z $$  (2.1)

 $$\displaystyle {\sigma_{zx} \over G} = {\partial u_x \over \partial z} + \theta y $$  (2.2)

Now, substituting equation (1.1) into equation (2.1), and similarly substituting equation (1.2) into equation (2.2) the result is

 $$\displaystyle {-G \theta z \over G} = {\partial u_x \over \partial y} - \theta z $$

 $$\displaystyle {G \theta y \over G} = {\partial u_x \over \partial z} + \theta y $$

In both equations, each $$\displaystyle G $$ cancels out. Also, $$\displaystyle \theta x $$ and $$\displaystyle \theta y $$ cancel out when subtracted from each side of the respective equation, thus leaving

 $$\displaystyle {\partial u_x \over \partial y} = 0 $$

 $$\displaystyle {\partial u_x \over \partial z} = 0 $$

Thus, we can conclude that there is no warping since

 $$\displaystyle u_x (y,z) = 0 $$

Flexural Shear Flow in Thin-Walled Cross Sections

Unsymmetrical thin walled cross sections
There will be a comparison drawn between cross sections with symmetry and cross sections of solids without symmetrical cross sections.

This is the generalized formula for a unsymmetrical cross section. (Equation 5.1)

$$ \int_{A_i}^{}{\frac{d\sigma _{xx}}{dx} dA}=-q_A$$

Case 1: Symmetrical about the Y-axis
The Stress factor is given under symmetrical cross sectional conditions as follows:

$$\sigma _{xx}=\frac{M_yZ}{I_y} $$

The shear flow factor $$q\;$$ is:

$$q(s)=\frac{-V_zQ_y}{I_y}$$

$$Q_y=\int_{A_S}^{}{ZdA}=Z_cA_s$$

Case 2:Unsymmetrical about the Y-axis
The Stress factor is given under unsymmetrical cross sectional conditions as follows:

$$\sigma _{xx}=(K_yM_z-K_{yz}M_y)y+(K_yM_y-K_{yz}M_z)z \;$$

$$K_y=\frac{I_y}{D}$$

$$K_{yz}=\frac{I_yz}{D}$$

$$K_z=\frac{I_Z}{D}$$

The stress in the x-direction can be solved for in matrix format.

$$\sigma _{xx}=\begin{vmatrix} y & z \end{vmatrix}_{1x2}\begin{vmatrix} K_y & -K_{yz}\\ -K_{yz} & K_z \end{vmatrix}\begin{Bmatrix} M_z\\M_y

\end{Bmatrix}\Leftrightarrow \sigma _{xx}=\begin{vmatrix} y & z \end{vmatrix}_{1x2}\begin{vmatrix} K_yM_Z & -K_{yz}M_y\\ -K_{yz}M_z & K_zM_y \end{vmatrix}$$

$$\sigma _{xx}=\begin{vmatrix} z & y \end{vmatrix}_{1x2}\begin{vmatrix} K_y & -K_{yz}\\ -K_{yz} & K_z \end{vmatrix}\begin{Bmatrix} M_y\\M_z

\end{Bmatrix}$$ It is more convenient for convention issues to use this version of the formula for $$\sigma _{xx}$$

Setting the moment about the z-axis to zero we are particularizing the to cross section of a symmetrical case

$$M_z=0\;$$

$$I_{yz}=0\Rightarrow D=I_yI_z$$

D being previously stated, $$K_y$$ and $$K_z$$ are as follows: $$K_y=\frac{1}{I_z}$$

$$K_y=0\;$$

$$K_z=\frac{1}{I_y}$$

The Shear flow for the particularize case is:

$$q(s)=-(K_yV_y-K_{yz}V_z)Q_z-(K_zV_z-K_{zy}V_y)Q_y\;$$

 $$Q_Z=\int_{A_s}^{}{ydA}\; \; \; \; \;

Q_y=\int_{A_s}^{}{zdA} $$

= MATLAB =

Airfoil

Above is a graph of $$\sigma_{xx}$$ at all y positions in the middle wing section on the top skin. Due to the graph being concave up, it can be seen that the load is compressive on the top skin of the wing section.

Above is a graph of $$\sigma_{xx}$$ at all y positions in the middle wing section on the bottom skin. Due to the graph being concave down, it can be seen that the load is tensile on the top skin of the wing section.

The average compressive stress was found to be 144.2 MPa.

Simply-Supported Boundary Conditions
A simply-supported rectangular plate under an in-plane distributed compressive load is to be evaluated:

Buckling mode shapes can be determined by the following equation: $$\displaystyle \psi (x,y) = c_{mn}\sin\left(\frac{m \pi x}{a}\right)\sin\left(\frac{n \pi y}{b}\right)\, \ {\rm for} \ m,n = 1, 2, 3, \ldots$$ 1

Perspective views of the buckling shape of a simply-supported rectangular plate for various cases of $$\displaystyle\ m$$ and $$\displaystyle\ n$$

Now considering the factor with $$\displaystyle\ x$$ as the variable in the expression for $$\displaystyle\ \psi_{mn}(x,y)$$, the period, $$\displaystyle\ T$$, is calculated as follows:

$$ \sin ( \frac{m \pi (x + T)}{a} ) = \sin ( \frac{m \pi x}{a} )$$

$$\frac{m\pi (x+T)}{a}=\frac{m\pi (x)}{a}+2\pi n$$ where $$\displaystyle\ n=1,2,3,...$$

Evaluating at $$\displaystyle\ n=1$$

$$\frac{m\pi (x+T)}{a}=\frac{m\pi (x)}{a}+2\pi$$

$$ x+T = x + \frac{2\pi a}{m\pi}$$

$$ T = \frac{2a}{m}$$

The values of $$\displaystyle\ T$$ along with plots of the function are calculated and displayed below for values of $$\displaystyle\ m=1,...,5$$ to show that $$\displaystyle\ m$$ is indeed the number of half wave-lengths generated.

The critical buckling load $$\displaystyle\ (P_x)_{cr}$$ is calculated as follows:

$$\displaystyle (P_x)_{cr}=k_c\frac{\pi ^2 D}{b}$$ 1

where $$\displaystyle k_c(m, a/b):=\left(\frac{m b}{a}+\frac{a}{m b}\right)^2$$

and $$\displaystyle D:=\frac{E h^3}{12 (1 - \nu^2)}$$

A graph of $$\displaystyle\ k_c$$ versus aspect ratio, $$\displaystyle\ a/b$$, for various values of $$\displaystyle\ m$$ is shown below

Going one step further, the critical buckling stress can be calculated by

$$\displaystyle (\sigma_{xx})_{cr}=\frac{(P_x)_{cr}}{bh}=k_c\frac{\pi ^2 D}{b^2 h}$$ 1

Clamped Boundary Conditions
Next a rectangular plate with clamped boundary conditions is to be examined

The critcal buckling stress under an in-plane compressive load can be calculated as follows:

$$\displaystyle\ (\sigma_{xx})_{cr}=K\frac{E}{1-\nu ^{2}}(\frac{t}{b})^2$$

where $$\displaystyle K$$ is a function different from $$\displaystyle k_c$$, and $$\displaystyle t \equiv h$$ the plate thickness.

The buckling shape is now calculated using

$$\displaystyle\psi (x,y)=\frac{c}{4}\left(1 - \cos \frac{2 \pi x}{a}\right)\left(1 - \cos \frac{2 \pi y}{b}\right)$$

Conclusions
The figure below shows the critical buckling stress as a function of aspect ratio for the top panel (under compressive load from prescribed conditions) of the NACA 2415 Airfoil under both simply-supported and clamped boundary conditions. Also shown in green is the average compressive stress in the panel. Since stringers support the panel, it is neither exactly simply-supported or clamped; therefore, the actual critical stress must fall within the envelope prescribed by the simply-supported and clamped cases. It can be seen that the average compressive stress of the panel does encounter the critical stress for the clamped case but not the simply-supported one; consequently, depending on how close the panel relates to the two cases, buckling is a possibility if the panel behaves closely to the clamped boundary condition.

= Team Contribution =

Christopher Pergola EAS4200C.F08.WIKI.A 19:27, 21 November 2008 (UTC)

Chris FontanaEas4200c.f08.WIKI.E 19:11, 17 November 2008 (UTC)

Jeff Huenink Eas4200c.f08.wiki.f 17:30, 18 November 2008 (UTC)

Braden Barnes Eas4200c.f08.wiki.b 19:51, 20 November 2008 (UTC)

Eric Viale Eas4200c.f08.wiki.c 04:07, 21 November 2008 (UTC)

Michael Rolle Eas4200c.f08.wiki.d 20:49, 21 November 2008 (UTC)

license: Creative Commons - Attribution Share-Alike - https://creativecommons.org/licenses/by-sa/4.0/

title: User:EAS4200C.F08.WIKI.A/Homework6


text: File:Bike p6.jpg

title: File:Bike p6.jpg


text:
User:Eml4500.f08.ateam.didomenico/hw 6

MATLAB: Eigenvalue Analysis of Rectangular Three-Bar and Four-Bar Truss Systems
This problem is a second attempt at plotting the eigenvectors for the truss problem described in Homework Report #4. The unstable three-bar truss system is depicted as the truss at the left of the of the following figure while the stable four-bar truss is shown at the right.

Unstable Three-Bar Truss System
The global stiffness matrix for the problem is assembled as follows.

Using the global stiffness matrix, the following MATLAB command can be used to determine the eigenvalues and vectors for the truss. These are displayed as follows with V corresponding to the eigenvectors and D corresponding to the eigenvalues.

Plotting Eigenvectors
By plotting the columns of the eigenvector matrix that correspond to the zero eigenvalues, a mode of the truss can be depicted. In this case, the zero eigenvalues fall in columns 1, 2, 3, 4, and 5 of the eigenvalue matrix, D. Therefore, the first five columns of the eigenvector matrix, V, are plotted similarly to plotting displacements. These values are simply added to the original nodal positions to give the mode shape.

The following code was appended to the original truss problem MATLAB file to plot the eigenvectors.

First Column Eigenvectors
This plot represents the eigenvectors corresponding to the first column of the eigenvector matrix, V.

Second Column Eigenvectors
This plot represents the eigenvectors corresponding to the second column of the eigenvector matrix, V.

Third Column Eigenvectors
This plot represents the eigenvectors corresponding to the third column of the eigenvector matrix, V.

Fourth Column Eigenvectors
This plot represents the eigenvectors corresponding to the fourth column of the eigenvector matrix, V.

Fifth Column Eigenvectors
This plot represents the eigenvectors corresponding to the fifth column of the eigenvector matrix, V.

Stable Four-Bar Truss System
The global stiffness matrix for the problem is assembled as follows.

Using the global stiffness matrix, the following MATLAB command can be used to determine the eigenvalues and vectors for the truss. These are displayed as follows with V corresponding to the eigenvectors and D corresponding to the eigenvalues.

Plotting Eigenvectors
By plotting the columns of the eigenvector matrix that correspond to the zero eigenvalues, a mode of the truss can be depicted. In this case, the zero eigenvalues fall in columns 1, 2, 3, and 4 of the eigenvalue matrix, D. Therefore, the first five columns of the eigenvector matrix, V, are plotted similarly to plotting displacements. These values are simply added to the original nodal positions to give the mode shape, just as in the previous truss system.

The following code was appended to the original truss problem MATLAB file to plot the eigenvectors.

First Column Eigenvectors
This plot represents the eigenvectors corresponding to the first column of the eigenvector matrix, V.

Second Column Eigenvectors
This plot represents the eigenvectors corresponding to the second column of the eigenvector matrix, V.

Third Column Eigenvectors
This plot represents the eigenvectors corresponding to the third column of the eigenvector matrix, V.

Fourth Column Eigenvectors
This plot represents the eigenvectors corresponding to the fourth column of the eigenvector matrix, V.

MATLAB: Electric Pylon Problem Analysis
The truss system in this problem is a electric-line pylon consisting of 91 elements. The MATLAB code for accurately assembling this truss system was obtained from Dr. Loc Vu-Quoc and can be found here or seen below.

The objective of this problem is to perform additional MATLAB programming to find specific properties of the truss system based on given specifications and an applied load P. P is applied vertically downward at the far right tip of the arm of the pylon. All elements are constructed of 300M steel and have identical cross-sectional areas. Properties for this material were given by Dr. Vu-Quoc and verified online here. They are displayed in the table below.

Running the provided code produces the following plot showing the structure of the electric pylon truss system.

Analysis Requirements

 * Scale the electric pylon so that its height is 60 m
 * Plot the undeformed shape (dotted line) and the deformed shape (solid line)
 * Compute the axial stress in each bar
 * Display the highest tensile stress and the highest compressive stress and the corresponding elements in which they occur
 * Indicate these elements with an arrow on the figure of undeformed and deformed shapes
 * Determine if the problem is statically determinant and provide justification through argument or MATLAB programming
 * Construct the lumped mass matrix of the electric pylon
 * Solve the generalized eigenvalue problem
 * Find the lowest three eigenpairs and, in 3 separate plots, plot the eigenvectors as deformed shapes (solid line) superposed onto the undedeformed shape (dotted line)
 * Find the 3 lowest vibrational periods of the electric pylon

MATLAB Analysis Code
The code below is the  file written to solve the problem as assigned. This program requires additional MATLAB functions to run. These functions are also detailed below.

Additional Functions
Similar to the five-bar truss system, the electric pylon truss system analysis MATLAB program requires additional functions to run. Some of hese functions are found at the course textbook website, Fundamental Finite Element Analysis and Applications. They include,  , and. Other essential functions include  and. These functions were independently written and cannot be found with the preceding functions. The details of these functions are outlined below.

This function generates the element stiffness matrix for each element. When the function is called it calculates the X and Y coordinate positions of the element ends from the coordinates passed in. It produces the element length from this information. Next it calculates the direction cosines of the element. Using the Young's modulus and cross-sectional area of the element, the element stiffness matrix, k, is assembled. The stiffness matrix is passed back to the.

This function solves the force-displacement relationship for each node, resulting in the displacements and reactions at each node. This information is passed back to the.

This function outputs the axial strain, axial stress, and axial force of each element.

This function returns the length of an element based on the element coordinates it receives. It was independently written and cannot be found with the preceding functions.

This function was written to output the reduced matrices, $$\bar{M}$$ and $$\bar{K}$$. It was independently written and cannot be found with the preceding functions.

Running Analysis Code
The code below is the complete results from running the MATLAB program. Each part is detailed in later sections of this homework report.

Plot of Undeformed and Deformed Truss System
The extended MATLAB code plots the deformation of the truss system as follows. To scale the electric pylon truss to be 60m tall, the nodal positions were multiplied by a factor of 8.0753701211306. It is important to note that the deformation has been magnified by a factor of 100 for viewing purposes.

Axial Stress of Elements
The axial stresses for each element are simply the second column of the results array. They are shown above as the result  when the MATLAB code is run. They are also shown here.

Highest Tensile Stress
As shown on the plot and results from running, the greatest tensile stress occurs in Element 81. The value of the stress is 9.0511x106 Pa.

Highest Compressive Stress
As shown on the plot and results from running, the greatest compressive stress occurs in Element 55. The value of the stress is 8.6957x106 Pa. It is important to note that the negative sign has been removed from this stress because it has been denoted as compressive stress. In the stress table above, negative stresses denote compressive stress. This aligns with common stress notation.

Lumped Mass Matrix
The lumped mass matrix was assembled as a diagonal matrix of the masses of both degrees of freedom for each node. To determine the masses, each elements length was multiplied by its cross-sectional area and its density. Once the global nodes of each end of the element were known, the half of the element's mass was assigned to both nodes. Continuing this for each element created the lumped mass matrix.

Eigenpair Analysis
The three smallest eigenvalues for the truss system were determined using the generalized eigenvalue problem shown below. In this relation $$\bar{K}$$ is the reduced stiffness matrix and $$\bar{M}$$ is the reduced lumped mass matrix.

$$\bar{K}\nu =\lambda \bar{M}\nu$$

By plotting the columns of the eigenvector matrix that correspond to the smallest eigenvalues, modes of the truss can be depicted. In this case, the smallest eigenvalues were 132.16, 2468.4 and 2940.3. Therefore, the columns of the eigenvector matrix, V, that correspond to these values are plotted similarly to plotting displacements. These values are simply added to the original nodal positions to give the mode shape.

The following code was appended to the original truss problem MATLAB file to plot the eigenvectors. It is important to note that the eigenvectors were multiplied by a factor of 10 for viewing purposes.

Smallest Eigenpair Plot
This plot represents the eigenvectors corresponding to the smallest eigenvalue.

Second Smallest Eigenpair Plot
This plot represents the eigenvectors corresponding to the second smallest eigenvalue.

Third Smallest Eigenpair Plot
This plot represents the eigenvectors corresponding to the third smallest eigenvalue.

Lowest Vibrational Periods of the Electric Pylon
To find the three lowest vibration periods of the electric pylon, the eigenvalues are set equal to the angular frequencies as shown in the equation below.

$$\lambda = \omega^2$$

The angular frequencies are converted to vibrational frequencies as follows.

$$f = \frac{\omega}{2\cdot \pi}$$

From the vibrational frequencies, the vibrational periods can be determined using the following equation.

$$T = \frac{1}{f} = \frac{2\cdot \pi}{\omega}$$

Following these relations reveals the three lowest vibrational periods, shown in the table below.

Statically Determinant?
The electric pylon structure is not statically determinate. While two nodes of the system are fixed, these two nodes act as the termination point of four elements. Thus, there are a total of four unknown reaction forces in a system with only three equations with which to attempt to solve them. This cannot be solved using the principles of statics.

It is important to note that, while similar in number of fixed nodes to previous systems analyzed, this system differs in that the fixed nodes contain more than one element.

title: User:Eml4500.f08.ateam.didomenico/hw 6


title: File:Fighw6d.jpg


text: User:Wikisaver62 (IMP)

title: User:Wikisaver62 (IMP)


text:
User:Eml4500.f08.RAMROD.F/HW6

= Mtg. 33 Wed November 12 =

Element stiffness matrix for element i
The Element stiffness matrix originates from the following equation:

$$\int_{x_{x}}^{x_{i+1}}{\left[N^1_iW_i+N^1_{i+1}W_{i+1} \right](EA)}\left[N^1_id_i+N^1_{i+1}d_{i+1} \right]dx$$

$$N^1_i=\frac{dN_i(x)}{dx}$$ Likewise for N_{i+1}

Note: $$u(x)=\begin{vmatrix} N_i(x) & N_{i+1}(x) \end{vmatrix}\begin{Bmatrix} d_i\\ d_{i+1} \end{Bmatrix}$$ and $$\frac{du(x)}{dx}=\begin{vmatrix} \dot{N}_i(x) & \dot{N}_{i+1}(x) \end{vmatrix}\begin{Bmatrix} d_i\\ d_{i+1} \end{Bmatrix}$$ where $$\mathbf{N}(x)_{1x2}=\begin{vmatrix} N_i(x) & N_{i+1}(x) \end{vmatrix}$$ and $$\mathbf{B}(x)_{1x2}=\begin{vmatrix} \dot{N}_i(x) & \dot{N}_{i+1}(x) \end{vmatrix}$$

Similiarly, the matrix $$\mathbf{W}(x)$$ can be written as:

$$\mathbf{W}(x)= \mathbf{N}(x)\begin{Bmatrix} W_i\\ W_{i+1} \end{Bmatrix}$$

$$\frac{d\mathbf{W}(x)}{dx}= \mathbf{B}(x)\begin{Bmatrix} W_i\\ W_{i+1} \end{Bmatrix}$$

Recalling the Degree of Freedom to establish these two in a matrix conventions:

$$ \mathbf{d}^{(i)}=\begin{Bmatrix} d_i\\ d_{i+1} \end{Bmatrix}\begin{Bmatrix} d^{(1)}_1\\ d^{(1)}_{2} \end{Bmatrix}$$ and $$ \mathbf{W}^{(i)}=\begin{Bmatrix} W_i\\ W_{i+1} \end{Bmatrix}\begin{Bmatrix} W^{(1)}_1\\ W^{(1)}_{2} \end{Bmatrix}$$

$$B=\int_{x_i}^{x_{i+1}}{\left[ (BW^{(i)})(EA)(Bd^{(i)})\right]dx}$$ $$B=W^{(i)}\cdot (K^{(i)}d^{(i)})$$

Using the transpose of matrix B to prove that the global stiffness matrix K will be listed as follows:

$$\mathbf{K}^{(i)}_{2x2}=\int_{x_i}^{x_{i+1}}{\left[ (B^{\tau } (x)(EA)(B(x)\right]dx}$$

which leads to the formation of the matrix K when EA is a constant being:

$$\mathbf{K}^{(i)}_{2x2}=\frac{EA}{L^{(i)}}\begin{bmatrix} 1 &-1 \\ -1 & 1 \end{bmatrix}$$

Transformation of Variable (coordinate) from $$x$$ to $$\tilde{x}$$
$$\tilde{x}=x-x_i$$ and $$d\tilde{x}=dx$$

so the B matrix in the new coordinate system can be expressed this way:

$$ \mathbf{K}^{(i)}_{2x2}=\int_{\tilde{x}=0}^{\tilde{x}=L^{(i)}}(B^{\tau } (\tilde{x})(EA){B(\tilde{x})d\tilde{x}}$$

which leads to the formation of the matrix K when EA and L is the length of element i is a constant being:

$$\mathbf{K}^{(i)}_{2x2}=\frac{EA}{L^{(i)}}\begin{bmatrix} 1 &-1 \\ -1 & 1 \end{bmatrix}$$

Here the area matrix and the elasticity matrix are:

$$A(\tilde{x})=N^{(i)}_1(\tilde{x})A_1+N^{(i)}_2(\tilde{x})A_2 $$

$$E(\tilde{x})=N^{(i)}_1(\tilde{x})E_1+N^{(i)}_2(\tilde{x})E_2$$ Giving us the Stiffness matrix K equaling: $$\mathbf{K}^{(i)}(\tilde{x})=\frac{\left[N^{(i)}_1(\tilde{x})E_1+N^{(i)}_2(\tilde{x})E_2 \right]\left[N^{(i)}_1(\tilde{x})A_1+N^{(i)}_2(\tilde{x})A_2 \right]}{L^{(i)}}\begin{bmatrix} 1 &-1 \\ -1& 1 \end{bmatrix}$$

title: User:Eml4500.f08.RAMROD.F/HW6


text:
User:Wikisaver62

Hi this is me. If contacting me then post a message on my wikipedia english userpage.

title: User:Wikisaver62


text: File:FEA Nov 7 Pic 6.jpg

title: File:FEA Nov 7 Pic 6.jpg


text:
User:Eas4200c.f08.radsam/HW6

 Comment: This page is identical to the previous one, except that the "Read and Report" section was added, which was accidentally omitted. Eas4200c.f08.radsam.a. 22:49, 6 December 2008 (UTC) =Homework #6=

=Class Notes=

November 3th - Plate Buckling
Plate buckling occurs when the applied axial load exceeds a critical value. Consider plate of dimensions:

a = dimension along the x-axis b = dimension along the y-axis

First, we look at a 1-D case:

$$[f] = \frac{F}{L}[A] = L^2\,\,$$

Therefore:

$$\left [ \frac{f}{A} \right ] = \frac{F}{L^3}$$ $$[ \sigma ] = \frac{F}{L^2} \Rightarrow \left [ \frac{\partial\sigma}{\partial x} \right ] = \frac{\frac{F}{L^2}}{L}\,$$ $$\epsilon = \frac{du}{dx}and\epsilon = \frac{\triangle L}{L}\Rightarrow [\epsilon ] = \frac{[du]}{[dx]} = \frac{L}{L} = 1\,$$

$$\nu = -\frac{\epsilon_{yy}}{\epsilon_{xx}}\Rightarrow [\nu] = 1\,$$

Eas4200c.f08.radsam.d 21:39, 21 November 2008 (UTC)

November 5 - Torsional Analysis
 Read and Report: Eqn 3.14

$$\frac{\partial \tau _{xz}}{\partial x} - \frac{\partial \tau _{yz}}{\partial y}=0$$  (Equation 3.14)

This equation shows that $$\tau _{yz}$$ and $$\tau _{xz}$$ are the only two nonvanishing stress components. Both are independent of z and are a reduction of the equations of equilibrium.

Looking back at the roadmap from earlier in the semester, recall the D-prandtl stress function:

$$\sigma _{yx} = \frac{\partial \phi }{\partial z}$$, $$\sigma _{zx} = -\frac{\partial \phi }{\partial y}$$

Notice that phi plays the role of a potential function. Thus, $$\sigma _{yx}$$ and $$\sigma _{zx}$$ are components of the "gradient" of phi with respect to y and z.

Also recall the scalar function f(x,y,z):

grad(f(x,y,z)) = $$\frac{\partial f}{\partial x}\hat{i}+\frac{\partial f}{\partial y}\hat{j}+\frac{\partial f}{\partial z}\hat{k}$$    where f is the potential function

By plugging this into our equations of equilibrium, we see that the equation is satisfied:

$$\frac{\partial }{\partial x}(\frac{\partial \phi }{\partial z})+\frac{\partial }{\partial x}(-\frac{\partial \phi }{\partial y}) = \frac{\partial^{2} \phi }{\partial y \partial z} - \frac{\partial^{2} \phi }{\partial z \partial y} = 0$$

This is true because phi is continuous and smooth, so the second derivative is interchangable. Thus, they are equal and subtracting them yields zero.

 Read and Report: Eqn 3.15 - 3.19

The same procedure can be followed for strain-displacement relations to put in terms of the Prandtl stress function:

$$\gamma _{xz} = \frac{\partial w}{\partial x}+\frac{\partial u}{\partial z}$$ and $$\gamma _{yz} = \frac{\partial w}{\partial y}+\frac{\partial v}{\partial z}$$

Thus, we obtain:

$$\gamma _{xz} = \frac{\partial w}{\partial x}-\theta y$$ and $$\gamma _{yz} = \frac{\partial w}{\partial y}+\theta x$$

We can now easily derive the following:

$$\frac{\partial \gamma _{yz}}{\partial x}-\frac{\partial \gamma _{xz}}{\partial y} = 2\theta$$

This is called the compatibility equation. Recalling what we proved earlier for the stress-strain relations:

$$\frac{\partial \tau _{yz}}{\partial x}-\frac{\partial \tau _{xz}}{\partial y} = 2G\theta$$

Finally, in terms of the Prandtl stress function, we end up with:

$$\frac{\partial^{2} \phi }{\partial y^{2}} + \frac{\partial^{2} \phi }{\partial z^{2}} = -2G\theta $$ where the left hand side of the equation is the LAPLACIAN of phi.

 Read and Report: Through end of Section 3.2

The following images show how the traction force and normal force vectors are defined.

 * The first image refers to a general 2D case (Figure 2)

 * The second image is an example of the traction force applied to an aircraft (Figure 3)

The vector $$t$$ vanishes because no loads are applied on the lateral surface: $${t} = [\sigma] {n}$$

this way, the stress vector may be evaluated on the lateral surface, knowing that $$n_z = 0$$, thus: $$\begin{Bmatrix} t_x\\ t_y\\ t_z \end{Bmatrix} = \begin{bmatrix} 0 &0 & \tau _{xz}\\ 0 &0 &\tau _{yz} \\ \tau _{xz} &\tau_{yz} &0 \end{bmatrix} \begin{Bmatrix} n_x\\ n_y\\ 0 \end{Bmatrix}\, $$

Therefore, we have: $$t_x = 0t_y = 0$$ $$t_z = \tau_{xz} n_x + \tau_{yz} n_y = \frac{\partial \phi }{\partial y}n_x - \frac{\partial \phi}{\partial y}n_y$$

If we look at Figure 4, we may easily derive: $$n_x = \sin \eta = \frac{dy}{ds}\,$$ $$n_y = \cos \eta = \frac{dx}{ds}$$

So we may express $$t_z\,$$ as:

$$t_z = \frac{\partial \phi}{\partial y}\frac{dy}{ds} + \frac{\partial \phi}{\partial x}\frac{dx}{ds} = \frac{d\phi}{ds}\,$$ So, the free boundary condition $$t_z = 0\,$$ is given by: $$\frac{d\phi}{ds} = 0 or  \phi = constant$$ on the lateral surface.

Note: For solid sections with a single contour boundary, this constant may be approximated to zero.

We are interested in the shear stresses $$\tau_{xz}~and~\tau_{yz}\,$$ and the resultant torque. Considering a small area $$dA = dx dy\,$$, the torque is:

$$dT = x \tau_{yz} dA - y \tau_{xz} dA\,$$ $$dT = \begin{pmatrix} -x\frac{\partial \phi}{\partial x} -y\frac{\partial \phi}{\partial y} \end{pmatrix} dA\,$$

Taking the integral:

$$T = -\int \int_{A}^{}{}(x\frac{\partial \phi }{\partial x}+y\frac{\partial \phi }{\partial y})dxdy$$

$$ = 2\int \int_{A}^{}{\phi dxdy}-\int \left[x\phi  \right]_{x1}^{x2}dy-\int \left[y\phi  \right]_{y1}^{y2}dx$$

Recall, phi vanishes on the boundary contour, so the last two terms in the equation vanish. Thus:

$$ = 2\int \int_{A}^{}{\phi dxdy} $$

Using our equations from before, combined with the warping equations developed in previous homework presenations, we can use the integration to find the torsion constant J:

$$J=\frac{T}{G\theta }$$

Eas4200c.f08.radsam.s 20:54, 21 November 2008 (UTC) Eas4200c.f08.radsam.d 21:40, 21 November 2008 (UTC)

November 10th - Torsional Analysis (Cont')
Homework 6: Find relationship between K and K_c. Verify K for a/b=1, m=1.

Torsional Analysis Continuation
 Case 1 : Thin-walled cross-section (closed) e.g. NACA airfoil.

$$\displaystyle \phi~=~C_o~on~S_o$$    (3.40)   (Fig 3.11)

$$\displaystyle \phi~=~C_1~on~S_1$$ (3.41)

 Case 2 :

$$\displaystyle \phi~=~0~on~S_o$$    (3.246)

Uniform bar with solid, circular cross-section (sec 3.3)

Look at page 74

Expression for T & J in terms of $$\displaystyle \phi$$

$$\displaystyle \phi(y,z)~=~C(\frac{y^2}{a^2}~+~\frac{z^2}{b^2}~-~1)$$

$$\displaystyle b~=~a$$(circle) $$\displaystyle \phi~=~0~on~S_o$$ $$\displaystyle \delta^2 \phi = -2G\Theta$$ $$C~=~-~\frac{1}{2}a^2G\Theta$$

November 12th

If we were to look at a cross section of a homogeneous circular bar, and place the origin of a coordinate system in its center...

$$x^2+y^2=a^2\,$$

Where a is the radius of the circular boundary.

Now, if we were to study the lateral surface of the bar, we would observe that the stress vector (traction) would be equal to zero. This means that the traction free boundary condition states that that the change of the stress function over the lateral surface is zero.

$$\frac{d\phi}{ds}=0\,$$

or

$$\phi=constant\,$$

This would mean that for single contour boundary-solid sections, the constant is arbitrary and thus can be zero. Hence the boundary condition can be defined on the lateral surface of the bar as:

$$\phi=0\,$$

Now that we have defined the boundary condition (stated above), we write the stress function in a form that satisfies this condition:

$$\phi = C(\frac{x^2}{a^2}+\frac{y^2}{a^2}-1)\,$$

Next we take the Prandtl stress function:

$$\frac{\partial ^2 \phi}{\partial x^2}+\frac{\partial ^2 \phi}{\partial y^2} = -2G\theta\,$$

Having obtained the compatibility equation above, we plug into it the stress function, which yields:

$$C=-\frac{1}{2}a^2G\theta\,$$

The stress function and the equation for C above, indicate the fact that the center of the circular cross section is the center of twist.

Recalling the definition of torque:

$$T= 2\iint_A \phi dx dy\,$$

We have the torque as:

$$T=2C\iint_A (\frac{x^2}{a^2} + \frac{y^2}{a^2}-1)dxdy\,$$

$$T=2C\iint_A (\frac{r^2}{a^2}-1)dA\,$$

$$T=2C(\frac{I_p}{a^2}-A)\,$$

Where:

$$I_p=\iint_A r^2dA= \frac{1}{2}\pi a^4\,$$

Where the above eq. corresponds to the polar moment of inertia and also:

$$A=\pi a^2\,$$ since that is the cross sectional are of the circular bar.

Since

$$a^2A=2I_p\,$$

Therefore:

$$T=-\frac{2CI_p}{a^2} = \theta GI_p\,$$

and

$$J = I_p\,$$

The shear stresses are:

$$\tau_{xz}=\frac{\partial \phi}{\partial y} = 2C\frac{y}{a^2}=-G\theta y\,$$

$$\tau_{yz}=-\frac{\partial \phi}{\partial x} = -2C\frac{x}{a^2} = G\theta x\,$$

It is also worth noting that the shear stress can be expressed as:

$$\tau=\frac{Tr}{J}\,$$

We can also show that for bars with circular cross-sections under torsion, no warping occurs:

$$\gamma_{xz}=\frac{\sigma_{xz}}{G}=\frac{\partial w_{x}}{\partial z} - \theta_{y}\,$$

$$\gamma_{yz}=\frac{\sigma_{yz}}{G}=\frac{\partial w_x}{\partial y} + \theta_{x}\,$$

Now recalling the shear stresses derived above:

$$\sigma_{xz}=-G\theta y\,$$

$$\sigma{yz}=G\theta x\,$$

We can see that the warping is zero:

$$w=0\,$$

Flexural Shear Flow in thin-walled sections

(Chapter 5) --> 5.1, 5.1.1, 5.1.2, 5.1.3, 5.3, 5.3.1, 5.3.2, 5.4

HW #5: Effects of M(x) on NACA airfoil

HW #6: Plate buckling, compressive in-plane load NACA airfoil

HW #7: Effects of V(x) on NACA airfoil, Plate Buckling, Shear Load.

Contribution by radsam.z

November 14th
 Flexural Shear Flow in Open, Asymmetric, Thin-Walled Sections

Sections 5.1, 5.1.1, 5.1.2

As $$\Delta x\,$$ approaches 0, the following equation can be used for either symmetric or asymmetric cross sections:

$$\int_{A_s}^{} \frac{d\sigma_{xx}}{dx} dA = -q_s$$

$$I_{yz}=0\,$$ due to cross section symmetry about either the y or z axis.

For symmetric cross sections, we have the following equations:

$$\sigma_{xx}= \frac{M_yz}{I_y}$$

$$q(s)=-\frac{V_zQ}{I_y}$$

$$Q=\int_{A_s}^{}z dA = z_c A_s$$

For asymmetric cross sections, we have:

$$\sigma_{xx} = (k_y M_z - k_{yz} M_y)y + (k_z M_y - k_{yz} M_z)z(5.3)\,$$

$$k_y= \frac {I_y}{I_yI_z-I_{yz}^2}\,$$

$$k_z = \frac {I_z}{I_yI_z-I_{yz}^2}\,$$

$$k_{yz}= \frac {I_{yz}}{I_yI_z-I_{yz}^2}\,$$

$$\sigma_{xx_{(1x1)}}=\begin{bmatrix}z & y\end{bmatrix}_{1x2}\begin{bmatrix}k_z & -k_{yz}\\ -k_{yz} & k_y\end{bmatrix}\begin{Bmatrix}M_y \\ M_z\end{Bmatrix}\,$$

Where the last two terms of the above product are: $$ \begin{Bmatrix} k_z M_y - k_{yz} M_z \\ -k_{yz} M_y + k_y M_z \end{Bmatrix}_{2x1} \,$$

We now particularize to a symmetric cross-section.

$$I_{yz}=0\,$$. [Why? See HW6 contribution below.] Consider $$M_z=0\,$$.

$$I_{yz}=0 => D=I_{y}I_{z}\,$$

$$k_y=\frac{1}{I_z}, k_{yz}=0, k_z=\frac{1}{I_y}\,$$

$$=>\sigma_{xx}=\frac{z M_{y}}{I_{y}}\,$$

For the asymmetric case:

$$q(s)=-(k_{y}V_{y}-k_{yz}V_{z})Q_{z}-(k_{z}V_{z}-k_{yz}V_{y})Q_{y}(5.5)\,$$

$$Q_{z}=\int_{A_{s}}ydA\,$$

$$Q_{y}=\int_{A_{s}}zdA\,$$

Stringer-web sections: Thickness of skin and spar webs is very small so neglect in components $$I_y\,$$, $$I_z\,$$, $$I_{yz}\,$$, $$Q_y\,$$, and $$Q_z\,$$. (Use only stringer areas.)

Symmetric with respect to y-axis:

Ex. 5.2, page 154 in book:

$$A_3=A_2\,$$

$$A_4=A_1\,$$

Asymmetric:

$$A_1=A\,$$

$$A_2=2A\,$$

$$A_3=3A\,$$

$$A_4=4A\,$$

Contribution by radsam.a

=Read and Report=

Chapter 3, Section 2 - Torsion of Uniform Bars
Torque is commonly encountered in an aircraft and is a very important stress component to consider. A torque is a applied about the longitudinal axis of a shaft or beam. A torque symbol (represented by a curved arrow) can contain multiple statically equivalent load distributions.

Two assumptions are made when deriving the equations for shaft analysis:
 * 1 - Cross sections remain planar upon the application of torque--thus, there is not out of plane deformation.
 * 2 - Diameters do not change upon deformation.

For non-circular cross sections, however, a different approach must be taken. There are two methods, this first being the Prandtl stress function method, and the second the Saint-Venant warping function method. Saint-Venant's warping function method will be derived as follows. Consider a straight shaft with a constant cross-section where equal and opposite torques, $$T\,$$, are applied at each end. The center of twist (COT) is defined as the location in which no in-plane displacements occur. The origin of the coordinate system is to be chosen at this location.

The twist angle per unit length at $$z\,$$ is given by:

$$\theta=\frac{\alpha}{z}\,$$

Where $$\alpha\,$$ is the twist or rotation angle at $$z\,$$ relative to the end at of the beam, at $$z=0\,$$.

It can be assumed that cross sections warp, but the cross section projections rotate on the x-y plane in a rigid manner. We now look at the point, $$P\,$$, on the edge of the cross-section located at $$z\,$$. Given a torque, $$P\,$$ moves to $$P'\,$$ over the angle, $$\alpha\,$$, which is considered to be small. Assuming that the cross-section at the base of the shaft doesn't move, the displacement components at $$P\,$$ are given by:

$$u=-r\alpha sin\beta=-\alpha y=-\theta zy\;$$

$$v=r\alpha cos\beta=\alpha x=\theta zx\;$$

where $$r\,$$ is the distance from the origin to point $$P\,$$.

The displacement in the $$z\,$$-direction, $$w\,$$, can be written:

$$w(x,y)=\theta \psi (x,y)\;$$

Which is the warping function. It should be noted that $$w\,$$ and $$\theta\,$$ do not depend on $$z\,$$. The given equations yield the following:

$$\epsilon_{xx}=\epsilon_{yy}=\epsilon_{zz}=\gamma_{xy}=0\,$$

Thus:

$$\sigma_{xx}=\sigma_{yy}=\sigma_{zz}=\tau_{xy}=0\,$$

So $$\tau_{yz}\,$$ and $$\tau_{xz}\,$$ are the two remaining stress components.

The equations of equilibrium can now be reduced to:

$$\frac{\delta \tau_{xz}}{\delta x}+\frac{\delta \tau_{yz}}{\delta y}=0$$

From Prandtl's stress function, $$\phi(x,y)\,$$, we have:

$$\tau_{xz}=\frac{\delta\phi}{\delta y}\,$$

$$\tau_{xz}=-\frac{\delta\phi}{\delta x}\,$$

In which both stress components automatically satisfy the equilibrium equations.

We then gather the following from the strain-displacement relations:

$$\gamma_{xz}=\frac{\delta w}{\delta x}+\frac{\delta u}{\delta z}\,$$

$$\gamma_{yz}=\frac{\delta w}{\delta y}+\frac{\delta v}{\delta z}\,$$

Which can be written as:

$$\gamma_{xz}=\frac{\delta w}{\delta x}-\theta y\,$$

$$\gamma_{yz}=\frac{\delta w}{\delta y}+\theta x\,$$

So the following torsion compatibility equation can be derived:

$$\frac{\delta \gamma_{yz}}{\delta x}-\frac{\delta \gamma_{xz}}{\delta y}=2\theta\,$$

This expressed in terms of the Prandtl stress function gives:

$$\frac{\delta^2 \phi}{\delta x^2}+\frac{\delta^2 \phi}{\delta y^2}=-2G\theta$$

Because of no loading on the lateral bar surface, the stress vector, $$t\,$$, can be evaluated as:

$$t_x=0\,$$

$$t_y=0\,$$

$$t_z=\tau_{xz}n_x+\tau_{yz}n_y=\frac{\delta \phi}{\delta y} n_x-\frac{\delta \phi}{\delta x} n_y$$

Where:

$$n_x=\frac{dy}{ds}\,$$ $$n_y=-\frac{dx}{ds}\,$$

Substituting into the stress vector equations gives:

$$t_z=\frac{\delta \phi}{\delta y}\frac{dy}{ds}+\frac{\delta \phi}{\delta x}\frac{dx}{ds}=\frac{d\phi}{ds}$$

$$t_z=0\,$$ is now expressed as:

$$\frac{d \phi}{ds}=0\mbox{or}\phi=constant$$

On the lateral surface. Choosing $$\phi=0\,$$ is allowable in that we are analyzing solid sections with a single contour boundary.

By integrating dT across the entire cross-section, the total resultant torque is given by:

$$T=-\int\int_A(x\frac{\delta \phi}{\delta x}+y\frac{\delta \phi}{\delta y})dxdy=2\int\int_A\phi dxdy-\int[x\phi]_{x_1}^{x_2}dy-\int[y\phi]_{y_1}^{y_2}dx$$

So:

$$T=2\int\int_A\phi dxdy$$

Finally, from:

$$J=\frac{T}{G\theta}\,$$

We get:

$$J=-\frac{4}{\nabla^2 \phi}\int\int\phi dxdy$$

Contribution by radsam.a

Chapter 3, Section 3 - Bars with Circular Cross-Sections
Looking at a uniform cross-section, the contour can be expressed as:

$$x^2+y^2=a^2\,$$

Where $$a\,$$ is the contour radius.

Then the stress function is then given by:

$$\phi=C(\frac{x^2}{a^2}+\frac{y^2}{a^2}-1)$$

Substituting into the compatibility equations results in:

$$C=-\frac{1}{2}a^2G\theta$$

So the center of the circular cross-section is also the center of twist!

The torque is given by:

$$T=2C(\frac{I_p}{a^2}-A)\,$$

Where:

$$I_p=\int\int_A r^2dA=\frac{1}{2}\pi a\,$$

gives the cross-section polar moment of inertia.

Also:

$$a^2 A=2I_p\;$$

So, substituting this into the equation for torque gives:

$$T=-\frac{2CI_p}{a^2}=\theta G I_p$$

And, $$J = I_p\,$$.

The shear stresses then given by:

$$\tau_{xz}=\frac{\delta \phi}{\delta y} = 2C\frac{y}{a^2}=-G\theta y$$

$$\tau_{yz}=-\frac{\delta \phi}{\delta x}=-2C\frac{x}{a^2}=G\theta x$$

And traction gives:

$$t_x=t_y=0\,$$

$$t_z=\tau_{xz}n_x+\tau_{yz}n_y\,$$

While:

$$n_x=cos\beta=\frac{x}{r}\,$$

$$n_y=sin\beta=\frac{y}{r}\,$$

So we end up with:

$$t_z=-G\theta\frac{xy}{r}+G\theta\frac{xy}{r}=0$$

Also:

$$n_z=0\,$$

And:

$$t_z=-G\theta{r}\,$$

So a little manipulation gives:

$$\tau=-t_z=G\theta r=\frac{Tr}{J}\,$$

We result in finding out that:

$$w=0\;$$

And we have confirmed that there is, indeed, no warping for circular cross-sections.

=Homework Problems and Solutions=

=Homework #5 Contributors= Eas4200c.f08.radsam.s 20:49, 21 November 2008 (UTC) Eas4200c.f08.radsam.d 19:42, 21 November 2008 (UTC) Eas4200c.f08.radsam.a. 21:20, 21 November 2008 (UTC) Eas4200c.f08.radsam.z 20:31, 20 November 2008 (UTC) Eas4200c.f08.radsam.r 14:09, 20 November 2008 (UTC) Eas4200c.f08.radsam.m 21:09, 21 November 2008 (UTC)

title: User:Eas4200c.f08.radsam/HW6


text:
A course in troll sockpuppets/Interview with a female troll

In response to criticisms that this course is male dominated, we located a female troll who agreed to be interviewed. Following is a transcript of the chat.

 * Durova: Hi, thanks for coming.
 * Krystal: Thank you so much for inviting me. It's wonderful to be here. [[Image:Love Heart SVG.svg|13px]]
 * Durova: So female trolls do exist?
 * Krystal: Oh yes. I'll send you my picture.
 * Durova: And you're typical?
 * Krystal: Well, I'm female.
 * Durova: You don't look like a troll in this picture.
 * Krystal: That's because I'm wearing my new Wonderbra. Do you think it fits?
 * Durova: Very...flattering...
 * Krystal: I'm lonely. Are you lonely too?
 * Durova: But I'm female.
 * Krystal: Oh that's great! I prefer girls.
 * Durova: This is just an interview. So how do you troll?
 * Krystal: No, I'm not a troll at all.
 * Durova: But you are a sockpuppet?
 * Krystal: This is my main account. I posed for a nude photo shoot last week and can't figure out which are the best pictures.  Would you help choose?
 * Durova: So if you're not a troll and you're not a sockpuppet, why did you agree to this interview?
 * Krystal: I need someone to talk to. Would you send me your picture?
 * Durova: Let's keep this on a professional level.

 * Troll: Hahahahahaha! You sick loser.  Hey, those pics I sent you had a virus.  My own creation.  You'll never get rid of it.

Instructor's note: it remains to be seen whether any actual female trolls exist.

Instructions
1. The basis for the female troll is a child's peach colored sock. Add a moderate amount of pillow filling to get a sense for the shape of the troll's head. Use that to size yarn for long blonde hair. 2. The female troll fits over the male troll's hand: one sockpuppet operates another. So fit the child's sock over the troll arm and adjust pillow stuffing to give the sock a basic puppet shape. 3. Cut ends of the hair and stitch into place with matching thread. To reinforce and prevent slippage, secure hair to scalp with a thin line of glue adjacent to the stitching. 4. Take a contrasting color sock the same size as the puppet sock and cut at the base of the ribbing. If desired, hem the cut end with elastic. 5. Fit the 'sweater' over the puppet, adding pillow stuffing for cleavage. 6. Create a pink female symbol and bend to puppet chest contours. Sew into place. 7. Sew buttons for eyes. Glue googly eyes over buttons. 8. Create a necklace of gold tone embroidery floss and size G glass beads. Three beads were used for this puppet, knotted to keep them at the front of the puppet. Attach a clasp if the puppet will be reused often. 9. Round off the head by folding down and stitching corners to head countour at back. 10. Cut pink felt to create lips. The lower lip on this puppet was glued into place. The upper lip was sewn.


text:
User:Eas4200c.f08.aero.e/HW6

=Week 11 - 11/3 through 11/7=

Plate Buckling
From Plate Buckling by Prof. Loc Vu-Quoc

This page talks about 4 different subjects. They are the compressive load and shear load, with simply-supported boundary conditions, and with Clamped boundary conditions.

Compressive Load
This image is of a rectangle of a and b under a compressive load. The thickness is represented by h.

#1 Simply-supported boundary conditions

The Buckling mode shapes is defined as $$\displaystyle \psi_{mn} (x,y)$$

$$ \displaystyle \psi (x,y) = c_{mn}\sin\left(\frac{m \pi x}{a}\right)\sin\left(\frac{n \pi y}{b}\right)\, \ {\rm for} \ m,n = 1, 2, 3, \ldots $$

where cmn is a constan, and (m,n) are the number of halfwave-lengths along (x,y), respectively.

The critical buckling load is defined as $$\displaystyle (P_x)_{cr}$$

$$ \displaystyle (P_x)_{cr} = k_c \frac{\pi ^2 D}{b} $$

Where kc and D, the plate bending stiffness, are

$$ \displaystyle k_c(m, a/b) = \left(\frac{m b}{a} + \frac{a}{m b}\right)^2 $$

$$ \displaystyle D = \frac{E h^3}{12 (1 - \nu^2)} $$

Then stated is the dimensions of;

[I] as length to the 4th [σ] as force over length squared, where є is 1 [EI] as force times length squared [D] as force time length, where $$\nu$$ is 1 [(Px)cr] as force, with kc is 1

Critical buckling stress is denoted below

$$ \displaystyle (\sigma_{xx})_{cr} = \frac{(P_x)_{cr}}{bh} = k_c\frac{\pi ^2 D}{b^2 h} $$

With the minimum, denoted with a star, critical load is

$$ \displaystyle (\sigma_{xx})_{cr}^\star = k_c^\star\frac{\pi ^2 D}{b^2 h} $$

where

$$  \displaystyle k_c^\star (a/b) =  \min_{m = 1, 2, \cdots} \ k_c (m, a/b) =  \min_{m = 1, 2, \cdots} \left(     \frac      {m b}      {a}      +      \frac      {a}      {m b}   \right)^2 =  \left(      \frac      {m^\star b}      {a}      +      \frac      {a}      {m^\star b}   \right)^2 $$

$$  \displaystyle m^\star (a/b) =  \underset{m=1,2, \cdots}{\rm argmin} \ k_c (m, a/b) =  \underset{m=1,2, \cdots}{\rm argmin} \left(     \frac      {m b}      {a}      +      \frac      {a}      {m b}   \right)^2 $$

And then for a simply-supported column, the critical load is

$$  \displaystyle (P_{cr})_{column} =  \frac {\pi^2 EI} {L^2} $$

#2 Clamped boundary conditions

An Excel file is given to find the critical stress given by the equation

$$  \displaystyle \sigma^\prime =  K   \frac{E}{1 - \nu^2} \left(     \frac{t}{b}   \right)^2 $$

from the book by Young & Budynas, Roark's formula for stress adn strain, McGraw-Hill, 7th ed., 2002, p.730.

In this equation t = h, and

$$  \displaystyle K  = k_c \frac{\pi^2}{12} $$

The estimated buckling load for a clamped plate close to being square from Timoshenko & Gere [1961], Theory of elastic stability, McGraw-Hill, p.386.

$$  \displaystyle \psi (x,y) =  \frac {c} {4}  \left(      1 - \cos \frac{2 \pi x}{a}   \right) \left(     1 - \cos \frac{2 \pi y}{b}   \right) $$

Shear Load
This image is the same rectangle under a shearing load.

#1 Simply-supported boundary conditions

The shear force per unit length is represented as

$$  \displaystyle N_{xy} =  \sigma_{xy} h $$

The critical buckling shear load is found by equating the strain energy and the work done by external shear load, and solving for the shear force.

The strain energy, U, and the work, W, are denoted as

$$  \displaystyle U  = \frac{D}{2} \frac{\pi^4 a b}{4} \sum_{m=1}^{\infty} \sum_{n=1}^{\infty} C_{mn}^2 \left(     \frac{m^2}{a^2}      +      \frac{n^2}{b^2}   \right)^2 $$

$$  \displaystyle W  = - 4 N_{xy} \sum_{m=1}^{\infty} \sum_{n=1}^{\infty} \sum_{p=1}^{\infty} \sum_{q=1}^{\infty} C_{mn} C_{pq} \frac {mnpq} {(m^2 - p^2)(n^2 - q^2)} $$

So

$$  \displaystyle (N_{xy})_{cr} =  -   \frac{ab D}{32} \frac {     \displaystyle \sum_{m=1}^{\infty} \sum_{n=1}^{\infty} C_{mn}^2 \left(        \frac{m^2 \pi^2}{a^2}	 +         \frac{n^2 \pi^2}{b^2}      \right)^2 }  {      \displaystyle \sum_{m=1}^{\infty} \sum_{n=1}^{\infty} \sum_{p=1}^{\infty} \sum_{q=1}^{\infty} C_{mn} C_{pq} \frac{mnpq}{(m^2 - p^2)(n^2 - q^2)} } $$

To find Cmn and Cpq to minimize the critical load, take the derivative with respect to Cmn and set it equal to zero.

Letting the aspect ratio of the plated be $$ \displaystyle \vartheta :=  \frac{a}{b} $$

$$  \displaystyle \lambda :=  -   \frac{\pi^2}{32 \vartheta} \frac{\pi^2 D}{b^2 h (\sigma_{xy})_{cr}} $$

The first 5 of these equations in matrix form are

$$  \displaystyle \left[ \begin{array}{lllll} \frac{\lambda (1 + \vartheta^2)^2}{\vartheta^2} &	 \frac{4}{9} &	 0	 &	 0	 &	 0	 \\	 \frac{4}{9} &	 \frac{16 \lambda (1 + \vartheta^2)^2}{\vartheta^2} &	 - \frac{4}{5} &	 - \frac{4}{5} &	 \frac{36}{25} \\	 0	 &	 - \frac{4}{5} &	 \frac{\lambda (1 + 9 \vartheta^2)^2}{\vartheta^2} &	 0	 &	 0	 \\	 0	 &	 - \frac{4}{5} &	 0	 &	 \frac{\lambda (9 + \vartheta^2)^2}{\vartheta^2} &	 0	 \\	 0	 &	 \frac{36}{25} &	 0	 &	 0	 &	 \frac{\lambda (9 + 9 \vartheta^2)^2}{\vartheta^2} \end{array} \right] \left\{ \begin{array}{l} C_{11} \\	 C_{22} \\	 C_{13} \\	 C_{31} \\	 C_{33} \end{array} \right\} =  \left\{ \begin{array}{l} 0	 \\	 0	 \\	 0	 \\	 0	 \\	 0     \end{array} \right\} $$

Witch is a symmetry diagonal matrix.

Case 1

To find the minimum critical buckling shear stress use the top left 2x2 matrix and set the determinant equal to zero. The eigenvalues are

$$  \displaystyle \lambda =  \pm \frac{1}{9} \frac{\vartheta^2}{(1 + \vartheta^2)^2} $$

Then the minimum critical buckling shear stress is

$$  \displaystyle (\sigma_{xy})_{cr}^{\star,2} =  \pm k_c^{\star,2} \frac{\pi^2 D}{b^2 h}  \, \ k_c^{\star,2} =  \frac{\pi^2}{32} \frac{9 (1 + \vartheta^2)^2}{\vartheta^3} $$

The 2 next to the star means 2x2 matrix and it is not dependent on the direction.

this will give a 15% error witch will increase with the increas in the aspect ratio.

Case 2

Using all 5 equations and setting the determinant to zero yields

$$  \displaystyle \lambda =  \left[ \frac{\vartheta^4}{81 (1 + \vartheta^2)^4} \left\{ 1	 +	 \frac{81}{625} +	 \frac{81}{25} \left(	   \frac{1 + \vartheta^2}{1 + 9 \vartheta^2}	 \right)^2 +	 \frac{81}{25} \left(	   \frac{1 + \vartheta^2}{9 + \vartheta^2}	 \right)^2 \right\} \right]^{1/2} $$

And

$$  \displaystyle (\sigma_{xy})_{cr}^{\star,5} =  \pm k_c^{\star,5} \frac{\pi^2 D}{b^2 h}  \, \ k_c^{\star,5} =  \frac{\pi^2}{32} \frac{1}{\lambda \vartheta} $$

#2 Clamped boundary conditions

The critical shear load is simply give by

$$  \displaystyle (\sigma_{xy})_{cr} =  K   \frac {E} {1 - \nu^2} \left(     \frac{h}{b}   \right)^2 $$

Where in Young & Budynas, Roark's formula for stress and strain, McGraw-Hill, 7th edition, 2002, p.732, a table of kc are given for aspect ratios from 1 to infinity. So K can be found.

It is also seen in experimental results a K = 4.1 for large aspect ratios. This may be due to imperfections in the experimental plate specimens.

--Eas4200c.f08.aero.E 06:38, 20 November 2008 (UTC)

Wednesday
$$\sigma_{xx}=\frac{d\phi}{dz}$$

$$\sigma_{zx} = -\frac{d\phi}{dy}$$

The above equations automatically satisfy the following equation: $$\frac{d\sigma_{yx}}{dy} + \frac{dz\sigma_{x}}{dz} = 0$$

In continuation from the previous lecture and lecture 27-2.

$$\sum_{3}^{i=1}{}\frac{2\sigma_{ij}}{dx_{i}} = 0$$ for j = 1, 2, 3

$$\begin{bmatrix} d\sigma_{ij}\\ dx_{i} \end{bmatrix} = \frac{F}{L^{3}} = \frac{Force}{Volume}$$

Recall roadmap on lecture 16-2. Prandtl stress funtion $$\phi$$

$$\sigma_{yx} = \frac{d\phi}{dz}$$

$$\sigma_{zx} = -\frac{d\phi}{dy}$$

$$\phi \Rightarrow $$Plays a role of a potential function.

$$(\sigma_{yx},\sigma_{zx})$$ comp. of "gradient" of $$\phi$$

Recall: Scalar funtion (x,y,z)

Gradient vector $$f(x,y,z) = \frac{df}{dx}i+\frac{df}{dy}j+\frac{df}{dz}k$$

f = potential

$$\frac{d}{dy}(\frac{d\phi}{dz})+\frac{d}{dz}(-\frac{d\phi}{dy})=\frac{d^{2}\phi}{dydz}-\frac{d^{2}\phi}{dzdy}=0$$

Since $$\phi$$ is continuous and smooth then the second dericative is interchangable

$$i.e. \Rightarrow\frac{d^{2}\phi}{dydz}=\frac{d^{2}\phi}{dzdy}$$

$$\frac{d^{2}\phi}{dy^{2}}+\frac{d^{2}\phi}{dz^{2}}=-2G\theta$$

Derive

$$[t]_{3x1} = [\sigma]_{3x3}{n}_{3x1}$$, equation 3.19 Chapter 2

{t} → Component of traction force

[σ] → Component of stress

{n} → Component of normal vector Ň

2-d Case on Page 29

Homework and reading assignments

1.) Read and report equation 3.14

It is known from a previous lecture that the stress-strain relations give the following:

$$\sigma_{xx}=\sigma_{yy}=\sigma_{zz}=\tau_{xy}=0$$

So, the only nonvanishing stress components are $$\tau_{yz}, and, \tau_{xz}$$ By thinking of the displacement field then it can be seen that these two stress components are independent of Z. If no body forces are considered then the equilibirum equations reduce to the following:

$$\frac{d\sigma_{yx}}{dy} + \frac{dz\sigma_{x}}{dz} = 0$$

2.) Read and report equation 3.15 to 3.19

Prandtl introduced the stress function Φ(X,Y), also τ_{zx} and τ_{yz} derived from Φ automatically satisfy the equations of equilibium. Now using the strain-displacements relations the following equation can be derived:

$$\frac{d\gamma_{yz}}{dx} - \frac{d\gamma_{xz}}{dy} = 2\theta$$

This is the compatibility equation for torsion, now thinking in terms of the prandtl stress-strain relations and the equation is derived

$$\frac{d^{2}\phi}{dx^{2}} + \frac{d^{2}\phi}{dy^{2}} = -2G\theta$$

2.) Read and report until the end of section 3.2

Section 3.2 mainly deals with torsion of uniform bars. Much like a bending moment a torque is a familiar load in          aircraft structures. Two assumptions are needed. One being that plane sections of a shaft remain plane and circular and the other that diameters remain straight after deformation are needed in deriving the deformation and stress fields.

There are two methods of sloving torsion of solid shafts, the Prandtl stress fuction method and the Saint-Venant warping function method. Both methods lead to the same results.

The torsion problem reduces the following two requirements; obtaining the stress fuction and satisfying the boundary conditions. Then the total torque is obtained by integrating dT over the entire cross-section.

--Eas4200c.f08.aero.arena 07:24, 20 November 2008 (UTC)

Friday
Example of $$\overrightarrow{t}$$ in aircraft

$$\displaystyle \overrightarrow{n} = n_y\overrightarrow{j} + n_z\overrightarrow{k}$$

$$\displaystyle ||n||= 1 $$

$$\overrightarrow{n} =$$ unit normal vector

$$ \displaystyle dz=ds cos \theta | n_y=cos\theta$$

$$\displaystyle dy=ds sin \theta | n_y=cos\theta$$

$$\displaystyle \sum F_y=0= -\sigma_{yy}*(dz*1)$$ (Where 1 is the unit depth along the x-dir)

 * $$\displaystyle -\sigma_{yz}*(dz*1) + t_y(ds*1)$$

$$\displaystyle 0 = -\sigma_{yy}*ds*n_y - \sigma_{yz}*ds*n_y + t_y*ds$$

Canceling ds we arrive at the following equation:

$$\displaystyle \sigma_{yy} n_y + \sigma_{yz}n_z = t_y $$ (1)

Note: $$[t_y] = \frac{F}{L^3}$$ => Units of force/area

$$\overrightarrow{t} =$$ traction vector (distributed surface force)

$$\displaystyle [t_y] = [\sigma] $$

$$\displaystyle t_z = \sigma_{yz}*n_y + \sigma_{zz}*n_z $$ (2)

Combining (1) and (2)

$$\displaystyle \begin{bmatrix} t_y \\ t_z \end{bmatrix}\, = \begin{bmatrix} \sigma_{yy} & \sigma_{yz} \\ \sigma_{zy} & \sigma_{zz} \end{bmatrix}\, \begin{bmatrix} n_y \\ n_z \end{bmatrix}\,$$ (3)

(3) can be generalized to the 3 dimensional case:

$$\begin{bmatrix} t_1 \\ t_2 \\ t_3 \end{bmatrix}\, = \begin{bmatrix} \sigma_{11} & \sigma_{12} & \sigma_{13} \\ \sigma_{21} & \sigma_{22} & \sigma_{23} \\ \sigma_{31} & \sigma_{32} & \sigma_{33} \end{bmatrix}\, \begin{bmatrix} n_1 \\ n_2 \\ n_3 \end{bmatrix}\,$$ (4)

$$\displaystyle t_i = \sum_{j=1}^{3} \sigma_{ij}n_{j} $$

--EAS4200C.F08.AERO.CLK04D 16:12, 21 November 2008 (UTC)

=Week 12 - 11/10 through 11/14=

Monday
For Homework Watch Ruth Goldberg Video Section 7.7 (Plate buckling)

Torsional Analysis Continued:

Case 1: Thin walled cross section (closed) e.g. NACA airfoil

Compare to fig 3.11

Case 2: Solid cross section

Now we consider the case of uniform bar with a solid, circular cross section. Refer to road map, expression for T and J in terms of $$\displaystyle \phi$$

$$ T = 2\int\int\limits_A\phi dA$$

$$J = \frac{-4}{\nabla^2\phi}\int\int\limits_A \phi dl $$

$$\phi(y,z) = C(\frac{y^2}{a^2} + \frac{z^2}{b^2} - 1)$$

(Note: b = 0 for a circle)

--EAS4200C.F08.AERO.CLK04D 16:10, 21 November 2008 (UTC)

Wednesday
The following graph will help us understand the derivation.

We know that:

$$r = \sqrt{x^{2}+y^{2}}$$

And, we also must recall that:

$$T = 2\oint_{A}^{}{}\phi dA=2C\left(\frac{J}{a^{2}}-A \right)$$ and $$J=\int_{A}^{}{}r^{2}dA = \frac{1}{2}\pi a^{4}$$

where $$A = \pi*a^{2}$$

so we must remember that: $$\begin{cases} \sigma _{yx}=\frac{\partial \phi }{\partial z}=-G\theta z \\ \sigma _{zx}=-\frac{\partial \phi }{\partial y}=G\theta y \end{cases}$$

-- Homework

Read and report through the end of section 3.3

 Report

Bars with circular cross-section

We consider a uniform bar with circular cross-section

We see that: $$x^{2}+y^{2}= a^{2}$$

And we assume the stress function of the following way: $$\phi = C\left(\frac{x^{2}}{a^{2}}+\frac{y^{2}}{a^{2}} -1\right)$$

Substituting the above expression into: $$\frac{\partial ^{2}\phi }{\partial x^{2}}+\frac{\partial ^{2}\phi }{\partial y^{2}}=-2G\theta $$

we get: $$C = -\frac{1}{2}a^{2}G\theta $$

Which solve the torsion problem and tells us that the center axis is the center of twist.

Now, we know that torque is: $$T=2\int_{}^{}{}\int_{A}^{}{}\phi dxdy$$

The using the above values and plug them in the Torque equation we get:

$$T = 2C\int_{}^{}{}\int_{A}^{}{}\left(\frac{x^{2}}{a^{2}}+\frac{y^{2}}{a^{2}}-1 \right)dxdy$$

That equals:

$$T = 2C\left(\frac{I_{p}}{a^{2}}-A \right)$$

where $$I_{p} = \int_{}^{}{}\int_{A}^{}{}r^{2}dA = \frac{1}{2}\pi a^{4}$$ and is the polar moment of inertia and A is the area of the cross-section.

And because $$A = \pi a^{2} $$ then $$t T = -\frac{2CI_{p}}{a^{2}} = \theta GI_{p}$$

So, the shear stresses are:

$$\begin{cases} \tau _{xz}=\frac{\partial \phi }{\partial y}= 2C\frac{y}{a^{2}} = -G\theta y \\ \tau _{yx}=-\frac{\partial \phi }{\partial x}=-2C\frac{x}{a^{2}}=G\theta x \end{cases}$$

Now, let's consider the figures below and concentrate on figure (a)

The stress vector are: $$t_{x}=t_{y}=0 $$ and $$  t_{z} = \tau _{xz}n_{x} + \tau _{yz}n_{y}$$

And since:

$$\begin{cases} n_{x}= cos\beta = \frac{x}{r}  \\ n_{y}=sin\beta = \frac{y}{r} \end{cases}$$

Then we get:

$$t_{z}=-G\theta \frac{xy}{r}+G\theta \frac{xy}{r}=0$$

We notice that the radial shear stress $$\tau _{rz}$$ vanishes because it is equal to $$t_{z}$$

Now, let's look at figure (b) of the above graph.

The unit vector on the surface is

$$n_{z}=0$$, $$n_{x}=sin\beta =\frac{y}{r}$$, $$ n_{y}=-cos\beta =-\frac{x}{r}$$

So, substituting the above equations as well as: $$\begin{cases} n_{x}= cos\beta = \frac{x}{r}  \\ n_{y}=sin\beta = \frac{y}{r} \end{cases}$$

into $$t_{z}=\tau _{xz}n_{x}+\tau _{yz}n_{y}=\frac{\partial \phi }{\partial y}n_{x}-\frac{\partial \phi }{\partial x}n_{y}$$

we get:

t_{z}=-G\theta r

We also see that on the z face the tangential shear stress is τz. And we end up getting:

$$\tau =-t_{z}=G\theta r$$

Where we use $$T=\theta =\theta GI_{p}$$ to get rid of θ

And finally, after algebra manipulation we get that $$w=u_{x}=0$$

Which tells us that for circular cross-sections bars there's no warping.

-- Homework

Using the equations: $$\begin{cases} \gamma _{yx} = \frac{\partial u_{x}}{\partial y}-\theta z \\ \gamma _{zx} = \frac{\partial u_{x}}{\partial z} + \theta y \end{cases}$$

and

$$\begin{cases} \tau _{yx} = \frac{\partial \theta }{\partial z} = \frac{2Cz}{a^{2}} = -G\theta z  \\ \tau _{zx} = -\frac{\partial \theta }{\partial y} =- \frac{2Cy}{a^{2}} = -G\theta y \end{cases}$$

And the δ-ε relations to show that

$$u_{x}\left(y,z \right) = 0$$

Solution

We know that:

$$\begin{cases} \tau _{yx} = \frac{\partial \theta }{\partial z} = \frac{2Cz}{a^{2}} = -G\theta z  \\ \tau _{zx} = -\frac{\partial \theta }{\partial y} =- \frac{2Cy}{a^{2}} = -G\theta y \end{cases}$$

Which can be rewritten as:

$$\begin{cases} \sigma_{yx} = -G\theta z&  \\ \sigma_{zx} = G\theta y& \end{cases} $$

Where we plug the σ values into the gamma values and we get:

$$\begin{cases} \gamma _{yx} = \frac{\partial u_{x}}{\partial y} - \theta z = \frac{\sigma _{yx}}{G} \Rightarrow -\theta z  \\ \gamma _{zx} = \frac{\partial u_{x}}{\partial z} + \theta y = \frac{\sigma _{zx}}{G} \Rightarrow \theta y \end{cases}$$

In a more simple way:

$$\begin{cases} \frac{\partial u_{x}}{\partial y} - \theta z =  -\theta z  \\ \frac{\partial u_{x}}{\partial z} + \theta y =  z\theta y \end{cases}$$

Solving for the partial part of the equation, we get:

$$\begin{cases} \frac{\partial u_{x}}{\partial y} = 0\\ \frac{\partial u_{x}}{\partial z} = 0 \end{cases} $$

Integrating the partial derivatives we get a constant and we realize that U is independent of x and y. -- Now, we proceed to start chapter 5

Motivation:

Now, if we draw a free body diagram of the above image in the following manner:

And the idea is to apply the principles developed in future clases for the NACA airfoil. Concepts and derivations for this graphs and chapter will be developed in the upcoming classes.

--Eas4200c.f08.aero.gabo28

Friday
$$\int_{A_s}\frac{d\sigma_{xx}}{dx}dA=-q_s$$(general for unsymmetrical cross section)

Symmetric

$$\sigma_{xx}=\frac{M_yz}{I_y}$$

$$q(s)=\frac{-V_zQ}{I_y}$$ $$Q=\int_{A_s}z dA=z_c A_s$$

Unsymmetrical Cross Section

$$\sigma_{xx}=(K_yM_-K_{yz}M_y)y+(K_zM_y-K_{yz}M_z)z$$

$$K_y=\frac{I_y}{D}$$ $$K_{yz}=\frac{I_{yz}}{D}$$ $$K_z=\frac{I_z}{D}$$

Where D is the determinant of I
 * $$\mathbf{\sigma_{xx}} = \begin{bmatrix}

y & z\end{bmatrix}$$ $$\begin{bmatrix} K_y & -K_{yz} \\ -K_{yz} & K_z\end{bmatrix}$$ $$\begin{bmatrix}

M_z\\ M_y\end{bmatrix}$$ Particularize to symmetric cross section, I_{yz}=0 consider $$M_z=0$$ $$I_{yz}=0$$ so $$D=I_yI_z$$ $$K_y=\frac{I_y}{D}=\frac{1}{I_z}$$, $$K_{yz}=0$$, $$K_z=\frac{I_c}{D}=\frac{1}{I_y}$$

so $$\sigma_{xx}=\frac{zM_y}{I_y}$$

Non-Symmetric $$q(s)= -(K_yV_y-K_{yz}V_z)Q_z -(K_zV_z-K_{yz}V_y)Q_y$$

stringer-web sections: thickness t of skin and spar webs are very small, so neglect in comp. $$I_y I_z I_{yz} Q_y Q_z$$ (use only areas of stringers)

EAS4200C.Fall08.AERO.Watlington.VG 21:57, 21 November 2008 (UTC)

=Matlab=

--Eas4200c.f08.aero.r 16:55, 21 November 2008 (UTC)--Eas4200c.f08.aero.r 16:55, 21 November 2008 (UTC)--Eas4200c.f08.aero.r 16:55, 21 November 2008 (UTC)--Eas4200c.f08.aero.r 16:55, 21 November 2008 (UTC)--Eas4200c.f08.aero.r 16:55, 21 November 2008 (UTC) --20:15, 21 November 2008 (UTC)20:15, 21 November 2008 (UTC)20:15, 21 November 2008 (UTC)~--Eas4200c.f08.aero.r 20:15, 21 November 2008 (UTC)

title: User:Eas4200c.f08.aero.e/HW6


text: Template:Using Wikiversity horizontal

title: Template:Using Wikiversity horizontal


text: File:Deform.JPG

title: File:Deform.JPG


text:
Template:Prevnext

A quick template for cycling through a sequence of pages:

Usage
E.G.

Where there is no previous (or no next) the parameter may be left empty. E.G.

Where the previous is the base page, use an "_". E.G.

Todo

 * 1) add pretty buttons
 * 2) add ifexist clauses
 * 3) add "up", first, last?
 * 4) If and when  attributes are enabled in MediaWiki, add   and   values, etc., to links.

title: Template:Prevnext


text: File:N.JPG

title: File:N.JPG


text:
Template:Chapter navigation

Related templates:
 * template:Chapter navigation top places small chapnav at top left
 * template:chapnavt shortcut to above
 * template:Chapter navigation places large chapnav bottom center
 * template:chapnav shortcut to above

This template adds:
 * 1) links to the previous and next pages of the book to the top of the page (position depends on stylesheet)
 * 2) links to the title page, previous page, and next page of the book to the bottom of the page
 * 3) the page to a category with the same name as the title page or optionally another category in its place, and sorts by the name of the subpage

Usage
Type at the bottom of each chapter in your book. If adding to the first chapter of a book, leave the first field blank. If adding to the last chapter, leave the second field blank. Pages will be added to a category of the book's title unless another category is specified in the third parameter.

can be used as a shortcut.

This can also be expanded to list the chapters of a book. See Learning Theories/Navigation and Chess/Navigation for examples.

Example
See Rhetoric and Composition/Narration for an example. The code  is used at the bottom of the page.

title: Template:Chapter navigation


text:
User:CGCgeo

This is the geography teacher. Culture  ES:cULTURA  PT:CULTURA

title: User:CGCgeo


text:
User:Eml4500.f08.ateam.nobrega/Homework Report Six

The Principle of Virtual Work: A Continuation of the Dynamics of an Elastic Bar
For the initial conditions, at t = 0:

 * $$ U(x,t=0) = \overline{U}(x) $$       The known functional displacement
 * $$ \dot{U}(x,t=0) = \overline{V}(x) $$       The known functional velocity

To create an equation for a multiple degree of freedom system, start with

 * $$ \frac{d}{dz}[EA \frac{du}{dx}] + f = m \ddot{u} $$    Equation (1)

This becomes,

 * $$ -kd + f = m \ddot{d} $$

or

 * $$ m \ddot{d} + kd = f $$    Equation (2)

To derive Equation 2 from 1 move all the terms to one side and take the integral with respect to x. This produces:

 * $$ \int_{x=0}^{x=L} w(x)(\frac{d}{dx}[EA \frac{du}{dx}] + f - m \ddot{u}) \,dx = 0 $$    Equation (3)

This is true for all w(x), the weighting function. It should be noted that producing Equation 3 from Equation 1 is a trivial solution, but producing Equation 1 from Equation 3 is not. Equation 3 can be re-written as,

 * $$ \int_0^L w(x)g(x) \,dx = 0 $$

Here g(x) is substituted into Equation 3 for

 * $$ \frac{d}{dx}[EA \frac{du}{dx}] + f - m \ddot{u} $$

Because Equation 3 holds true for all w(x), select w(x) such that it is equal to g(x). Doing so allows Equation 3 to become:

 * $$ \int_0^L g^2(x) \,dx = 0 $$

for $$ g(x) \ge 0 $$

Integration by Parts
When integrating an argument with two functions dependent upon the same variable, you must use the method of integration by parts.

 * $$ \int r(x)s(x)dx

(rs)'=r's + rs'

r'=\frac{dr}{dx}

s'=\frac{ds}{dx}

\int (rs)' = \int r's + \int rs'

\int (rs)' = rs - \int rs' $$

(Recall continuation of Principle of Virtual Work Equation 3 pg 29-1 also featured above)

First term: $$ r(x)=EA\frac{du}{dx}

s(x)=w(x) $$

by integration by parts:

 * $$ \int_{x=0}^{x=L} w(x)(\frac{d}{dx}[EA \frac{du}{dx}])dx = WEA\frac{du}{dx} - \int_{0}^{L} \frac{dw}{dx}(EA)[\frac{du}{dx}]dx

$$

The negative sign in front of the integral on the left hand of the equation corresponds to the negative sign in the equation:

$$ F-kd = M\ddot{d} $$

Elastic Beam Problem
Boundary Conditions: At x=0, select w(x) such that w(0)=0 (i.e. Kinematically Admissible)

The motivation is the discrete PVW applied to the equation below (also found on pg 10-1)

$$ W([k] \begin{bmatrix}     d_3 \\      d_4 \end{bmatrix} -F) = 0 $$

$$ F^{T}= \begin{bmatrix} F_{1} & F_{2} & F_{3} & F_{4} & F_{5} & F_{6} \end{bmatrix} $$

Where F3 and F4 are known and F1, F2, F5 and F6 are unknown reactions. Since W can be selected arbitrarily, select W such that W1=W2=W5=W6=0 so to eliminate equations involving unknown reactions. Thus, eliminating rows 1,2,5 and 6.

$$ kd=F $$

Where k is a 2x2 matrix, d is a 2x1 matrix and F is a 2x1 matrix. Note: W(kd-F)=0 for all values of W.

Unknown Reaction:

$$ N(0,t)=EA(0)\frac{du}{dx}(0,t) $$

Continue PVW with functions W(L) and F(t)

$$ -\int_{0}^{L}\frac{dW}{dx}(EA)\frac{du}{dx}dx + \int_{0}^{L} W(x)[f-m\ddot{u}]dx = 0 $$

For all W(x) such that W(0)=0

$$ \int_{0}^{L}w(m\ddot{u})dx + \int_{0}{L}\frac{dW}{dx}EA\frac{du}{dx}dx = W(L)F(t) + \int_{0}^{L}W f dx $$

For all W(x) such that W(0)=0.

Stiffness Terms:

Assume displacement u(x) for xi <= x <= xi+1
 * (i.e. x belongs to [xi, xi+1])

Motivation for linear interpolation of u(x): 2 Bar Truss

The deformed shape of the truss is a straight line as well because of the implicit assumption of linear interpolation of displacement between 2 nodes.

Consider the case where there are only axial displacements (i.e. 0 translational displacement).

Question: Express u(x) in terms of di=u(xi) and di+1=u(xi+1) as a linear function in x (i.e. linear interpolation)

u(x) = Ni(x)di + Ni+1(x)di+1

where Ni(x) and Ni+1(x) are linear functions of x.

N i+1 (x)= (x - xi) * (xi+1 - xi)-1

Honesty, Imagination, Ethics
The following quotes are excerpts from Dr. L. Vu-Quoc main homework page for the Finite Element Analysis: EML 4500 course at the University of Floida.

Continuation PVW to Discrete PVW
Langrangian Interpolation

Motivation for the form of $$ N_{i}(x) $$ and $$ N_{i+1}(x) $$

1. $$ N_{i}(x) $$ and $$ N_{i+1}(x) $$ are linear (straight lines). Thus any linear combination of $$ N_{i} $$ and $$ N_{i+1} $$ is also linear and in particular the expression for $$ u(x) $$ on p.31-3.

$$ N_{i}(x) = \alpha_{i} + \beta_{i}(x) $$, where $$ \alpha_{i} $$ and $$ \beta_{i} $$, are real numbers

$$ N_{i+1}(x) = \alpha_{i+1} + \beta_{i+1}(x) $$, where $$ \alpha_{i+1} $$ and $$ \beta_{i+1} $$, are real numbers

$$ N_{i}d_{i}+N_{i+1}d_{i+1} = (\alpha_{i}+\beta_{i}x)d_{i}+(\alpha_{i+1}+\beta_{i+1}x)d_{i+1} $$, is clearly a linear function in x

2. Recall the equation for $$ u(x) $$

$$ u(x_i) = N_{i}(x_i)d_i+N_{i+1}(x_i)d_{i+1} $$

Where $$ N_{i}(x_{i})d_i = 1 $$ and $$ N_{i+1}(x_{i})d_{i+1} = 0 $$

Continueing the FEM via PVW
From a previous lecture on meeting 31 about the linear interpolation of U(x), W(x) can be shown by applying the same principle.

i.e.

$$ W(x)=N_{i}(x)W_{i}+N_{i+1}(x)W_{i+1} $$

Now the element stiffness matrix for element I can be written as:

$$ \beta = \int_{x_{i}}^{x_{i+1}}{[N_{i}'W_{i}+N_{i+1}'W_{i+1}](EA)} $$

with

$$[N_{i}'W_{i}+N_{i+1}'W_{i+1}]=W'(x)$$

$$[N_i'd_i+N_{i+1}'d_{i+1}]=U'(x)$$

$$ N_i':=\frac{dN_i(x)}{dx} $$

Likewise for $$N_{i+1}'$$

NOTE

$$ u(x)=\begin{bmatrix} N_i(x) & N_{i+1}(x) \end{bmatrix} \begin{Bmatrix} d_i\\ d_{i+1} \end{Bmatrix} =\mathbf{N(x)_{(1x1)}} \begin{Bmatrix} d_i\\ d_{i+1} \end{Bmatrix} $$

$$ \frac{du(x)}{dx}=\begin{bmatrix} N_i'(x) & N_{i+1}'(x) \end{bmatrix} \begin{Bmatrix} d_i\\ d_{i+1} \end{Bmatrix} =\mathbf{B(x)} \begin{Bmatrix} d_i\\ d_{i+1} \end{Bmatrix} $$

Similarly W(x) can be shown as:

$$ W(x)=\mathbf{N(x)}\begin{Bmatrix} W_i\\ W_{i+1} \end{Bmatrix} $$

$$ \frac{dW(x)}{dx}=\mathbf{B(x)}\begin{Bmatrix} W_i\\ W_{i+1} \end{Bmatrix} $$

Now, let us recall the element degrees of freedom as discussed in a previous lecture

therefore

$$ \begin{Bmatrix} d_i\\ d_{i+1}

\end{Bmatrix} = \begin{Bmatrix} d_1^{(i)}\\ d_2^{(i)} \end{Bmatrix} = \mathbf{d^{(i)}} $$

$$ \begin{Bmatrix} W_i\\ W_{i+1}

\end{Bmatrix} = \begin{Bmatrix} W_1^{(i)}\\ W_2^{(i)} \end{Bmatrix} = \mathbf{W^{(i)}} $$

The element stiffness matrix can now be written as

$$ \beta =\int_{x_{i}}^{x_{i+1}}{{BW}^{(i)}_{(1x1)}(EA)_{(1x1)}{Bd}^{(i)}_{(1x1)}}dx = {W}^{(i)}\cdot({k}^{(i)}{d}^{(i)}) $$

$$ \beta =\int_{x_i}^{X_{i+1}}(EA)({BW^{(i)}})_{(1x1)}\cdot({Bd^{(i)}}){(1x1)}dx

$$

$$ ({BW^{(i)}})^T({Bd^{(i)}}) $$

where

$$ {BW^{(i)T}} = {W^{(i)T}B^T} = {W^{(i)} \cdot B^T} $$

$$ \beta= \mathbf{W^{(i)}}\cdot\int \mathbf{B^T}EA\mathbf{B}dx\mathbf{d^{(i)}} $$

$$ \mathbf{k^{(i)}}_{(2x2)}=\int_{x_i}^{x^{i+1}}{\mathbf{B(x)^T}}_{(2x1)}(EA)_{1x1}\mathbf{B(x)}_{(1x2)}dx $$

If we consider $$ EA $$ to constant the matrix $$ {k^{(i)}} $$ can be written as

$$\textbf{k}^{(i)} = \frac{EA}{L^{(i)}}\begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix}$$

The next step is the transformation of variable coordinate system from:

$${x}$$ to $$ \tilde{x} $$

Given

$$ d\tilde{x}=dx $$

$$ \tilde{x}:=x-x_i $$

Now $$ {k^{(i)}} $$ can be written as

$$ \mathbf{k^{(i)}}=\int_{\tilde{x}=0}^{\tilde{x}=L^{(i)}}{\mathbf{B^T(\tilde{x})}}EA(\tilde{x})\mathbf{B(\tilde{x}})d(\tilde{x}) $$

The following problem outlined by the image below can be solved using the transformed coordinate system.

Given:

With the following functions for the modulus elasticity and cross sectional area $$ {k^{(i)}} $$ matrix can be written as

$$A(\tilde{x}) = N_1^{(i)}(\tilde{x})A_1 + N_2^{(i)}(\tilde{x})A_2$$

$$E(\tilde{x}) = N_1^{(i)}(\tilde{x})E_1 + N_2^{(i)}(\tilde{x})E_2$$

$$ \mathbf{k^{(i)}}=\int_{\tilde{x}=0}^{\tilde{x}=L^{(i)}}{\mathbf{B^T(\tilde{x})}}{[N_1^{(i)}(\tilde{x})A_1 + N_2^{(i)}(\tilde{x})A_2]}{[N_1^{(i)}(\tilde{x})E_1 + N_2^{(i)}(\tilde{x})E_2]}(\tilde{x})\mathbf{B(\tilde{x}})d(\tilde{x}) $$

$$N_{1}^{(i)}(\tilde{x})=1-\frac{\tilde{x}}{L^{(i)}}= \begin{cases} 1 & \text{ at } \tilde{x}=0 \\ 0 & \text{ at } \tilde{x}=L^{(i)} \end{cases}$$

$$N_{2}^{(i)}(\tilde{x})=\frac{\tilde{x}}{L^{(i)}}= \begin{cases} 0 & \text{ at } \tilde{x}=0 \\ 1 & \text{ at } \tilde{x}=L^{(i)} \end{cases} $$

Homework:

Book p. 159 ; Set $$ E_{1}=E_{2}=E_{3}$$, let $$A(\tilde{x})$$ be linear as previously defined, and obtain $$\mathbf{k}^{(i)}$$ and compare to expression given in the book

Where:

$$A(\tilde{x})=N_{1}^{(i)}(\tilde{x})A_{1}+N_{2}^{(i)}(\tilde{x})A_{2}$$

and

$$E(\tilde{x})=N_{1}^{(i)}(\tilde{x})E_{1}+N_{2}^{(i)}(\tilde{x})E_{2}$$

Compare general $$\mathbf{k}^{(i)}$$ from before to $$\mathbf{k}_{avg}^{(i)}$$ where,

$$ \frac {(E_{1}+E_{2})}{2}\frac{(A_{1}+A_{2})}{2}\frac{1}{L^{(i)}} \begin{bmatrix} 1 &-1 \\ -1 & 1 \end{bmatrix}=\mathbf{k}_{avg}^{(i)}$$

Remember: Recall the mean value Theorem (MVT) and its relation to the centroid:

Mean Value Theorem: for $$\bar{x}\epsilon [a,b]$$,     $$a\leq \bar{x}\leq b$$,

$$ \int_{x=a}^{x=b}{f(x)dx}=f(\tilde{x})[b-a] $$

Where the centroid is defined as:

$$\int_{A}^{}{xdA} = \bar{x} \int_{A}^{}{dA} = \bar{x}A$$

Therefore:

$$\int_{x=a}^{x=b}{f(x)g(x)dx} = f(\bar{x})g(\bar{x})[b-a]$$ where  $$ a\leq \bar{x}\leq b$$

In general:

$$f(\bar{x})\neq \frac{1}{b-a} \int_{a}^{b}{f(x)dx}$$

$$g(\bar{x})\neq \frac{1}{b-a} \int_{a}^{b}{g(x)dx} $$

The Electric Pylon
Code was provided that creates a two dimensional truss system of an electric pylon. This code was to be scaled to 60 meters and loaded with a downward 1000 N force on the right most tip. The stresses in the element were to be calculated and the largest compressive and tensile stresses identified in the truss. The mass of each element needed to be taken into account with the mass of each element being redistributed to its end most points. The lowest eigenvectors were then plotted as well as the lowest periods of vibration.

Functions Needed to Load and Solve the Truss System
Length.m

The Length function calculates the element length based on the elements coordinates. function L = Length(coord) x_1 = coord(1,1); x_2 = coord(2,1); y_1 = coord(1,2); y_2 = coord(2,2); L = sqrt((x_2 - x_1)^2 + (y_2 - y_1)^2); %calculates element length using coordinates

ReduceM.m

This function takes the mass matrix and reduces it. function Mr = ReduceM(Mass, F, debc, ebcVals) %calculates the reduced matrix dof = length(F); df = setdiff(1:dof, debc); Mr = Mass(df, df);

ReduceS.m

This function reduces the stiffness matrix. function Kr = ReduceS(K, F, debc, ebcVals) %Reduces the stiffness matrix dof = length(F); df = setdiff(1:dof, debc); Kr = K(df,df);

PlaneTrussElement.m

This function calculates the length of the bar and the director cosines as well as the elemental stiffness matrix before returning it to the call function. function k = PlaneTrussElement(e, A, coord) % PlaneTrussElement(e, A, coord) % Generates stiffness matrix for a plane truss element % e = modulus of elasticity % A = area of cross-section % coord = coordinates at the element ends x1=coord(1,1); y1=coord(1,2); x2=coord(2,1); y2=coord(2,2); L=sqrt((x2-x1)^2+(y2-y1)^2); ls=(x2-x1)/L; ms=(y2-y1)/L; k = e*A/L*[ls^2, ls*ms,-ls^2,-ls*ms; ls*ms, ms^2,-ls*ms,-ms^2; -ls^2,-ls*ms,ls^2,ls*ms; -ls*ms,-ms^2,ls*ms,ms^2];

NodalSoln.m

This function uses the global stiffness matrix to calculate the number of global degrees of freedom. It calculates the global displacements and the reactions seen in the bars. function [d, reactions] = NodalSoln(K, R, deg, num) dof = length(R); diff = setdiff(1:dof, deg); Ksub = K(diff, diff); Rsub = R(diff) - K(diff, deg)*num; value = Ksub\Rsub; d = zeros(dof,1); d(deg) = num; d(diff) = value; reactions = K(deg,:)*d - R(deg);

PlaneTrussResults.m

This function creates the transform matrix T after calculating the length of the bars and the director cosines. This transform matrix turns the global displacements into axial displacements so it can compute the axial forces, stress, and strain acting within the bar. function results = PlaneTrussResults(e, A, coord, disps) % results = PlaneTrussResults(e, A, coord, disps) % Compute plane truss element results % e = modulus of elasticity % A = Area of cross-section % coord = coordinates at the element ends % disps = displacements at element ends % The output quantities are eps = axial strain % sigma = axial stress and force = axial force. x1=coord(1,1); y1=coord(1,2); x2=coord(2,1); y2=coord(2,2); L=sqrt((x2-x1)^2+(y2-y1)^2); ls=(x2-x1)/L; ms=(y2-y1)/L; T=[ls,ms,0,0; 0,0,ls,ms]; d = T*disps; eps= (d(2)-d(1))/L; sigma = e.*eps; force = sigma.*A; results=[eps, sigma, force];

Results
Highest Tensile Stress : 9.05 MPa

Element with Highest Tensile Stress: 81

Highest Compressive Stress: -8.70 MPa

Element with Highest Compressive Stress: 55

Lowest three vibrational periods :    0.5465 secs, 0.1265 secs , 0.1159 secs

Deformed Truss

Problem Statically Determinant?
This problem is not statically determinate. Each element has two degrees of freedom for each node, or a total of four degrees of freedom per element. With only three static equations, the problem is statically indeterminate. The two force member method does not apply as the forces cannot be assumed to act only in the direction of member.

Eigenvalue Plot Deformation

Eigenvalue Plot Deformation

Contributing Team Members
Eml4500.f08.ateam.nobrega 01:19, 21 November 2008 (UTC)

Eml4500.f08.ateam.boggs.t 01:37, 21 November 2008 (UTC)

Eml4500.f08.ateam.carr 02:50, 21 November 2008 (UTC)

Eml4500.f08.ateam.shah 18:21, 21 November 2008 (UTC)

Eml4500.f08.ateam.mcnally 18:51, 21 November 2008 (UTC)

Eml4500.f08.ateam.rivero 19:31, 21 November 2008 (UTC)


text:
User:Eml4500.f08.gravy.mmm/HW6

Rube-Goldberg device by (Walter Benjamin which is a German name so the W sounds like a V and the J sounds like a Y much like the Volkswagen where the V sounds like a F and the W sounds like a V) is a device which performs a simple task in a complicated manner. More information can be obtain here

"Honesty, imagination,and ethics" are important when it comes to writing these homework reports. Because everyone post their homework online it is important to remain honest and abide by your ethics and do your own work. When writing up the reports it takes imagination to use narrative in the reports.

From Continuous Principal of Virtual Work to Discrete Principal of Virtual Work
Discrete Principal of Virtual Work can be achieved using Lagrangian interpolation which is the motivation for the formation of Ni(x) and Ni+1(x).

1) Ni(x) and Ni+1(x) are linear straight lines, thus any linear combination of Ni and Ni+1 is also linear and in particular the expression for u(x) shown previously in this report.

$$N_i(x)= \alpha_i + \beta _ix$$, with ($$\alpha_i + \beta _i$$) are numbers

Linear combination of Ni and Ni+1
The Linear combination of Ni and Ni+1 is calculated as the following:

$$\ N_id_i + N_{i+1}d_{i+1}= (\alpha_i + \beta _ix)d_i +(\alpha_{i+1} + \beta _{i+1}x)d_{i+1} = (\alpha_{i}d_{i} +\alpha_{i+1}d_{i+1}) + (\beta _{i}d_{i}+\beta _{i+1}d_{i+1}) \,$$

which is clearly a linear function of x.

2)Recall equation for u(x)(the interpolation of u(x)) from previous in this report

$$\ u(x_i)= \underbrace{N_i(x_i)d_i}_1 + \underbrace{N_{i+1}(x_i)d_{i+1}}_0 = d_i \,$$

Prove of equation $$u(x_{i+1})= d_{i+1}$$
The same interpolation can be used for w(x) so     $$\ w(x)= N_i(x)w_i+ N_{i+1}(x)w_{i+1} \,$$

The Element Stiffness matrix for element i
The Element Stiffness matrix for element i can be developed using the steps and equations which follow.

$$ \beta= \int_{x_i}^{x_{i+1}} [ \underbrace{N'_iW_i + N'_{i+1}W_{i+1}}_{W'(x)} (EA) \underbrace{N'_id_i + N'_{i+1}d_{i+1}}_{u'(x)}]dx $$

Where $$ N^{\prime}_i:= \frac{d N_i(x)}{dx}$$ Likewise for $$N^{\prime}_{i+1}$$ and written in matrix form for simplification.

$$ u(x)= \underbrace{\begin{bmatrix} N_{i}(x) && N_{i+1}(x) \end{bmatrix}}_{\vec{N}_{1x2}(x)} \begin{Bmatrix} d_i \\ d_{i+1}  \end{Bmatrix}_{2x1}  $$

Recall the element degrees of freedom

$$ \begin{Bmatrix} d_{i} \\ d_{i+1} \end{Bmatrix}= \begin{Bmatrix} d^{(i)}_1 \\ d^{(i)}_{2}  \end{Bmatrix}= \mathbf{d}^{(i)}$$

$$ \begin{Bmatrix} W_i \\ W_{i+1} \end{Bmatrix}= \begin{Bmatrix} W^{(i)}_1 \\ W^{(i)}_{2}  \end{Bmatrix}= \mathbf{W}^{(i)}$$

$$ \beta= \int_{x_i}^{x_{i+1}} \underbrace{(\mathbf{BW}^{(i)}) }_{1x1} \underbrace{(EA)}_{1x1} \underbrace{(\mathbf{Bd}^{(i)}}_{1x1}dx = \mathbf{W}^{(i)} \cdot (\mathbf{k}^{(i)} \mathbf{d}^{(i)}) $$

$$ \beta= \int_{x_i}^{x_{i+1}}(EA) \underbrace{(\mathbf{BW}^{(i)}) }_{1x1}\cdot  \underbrace{(\mathbf{Bd}^{(i)})}_{1x1}dx = \int_{x_i}^{x_{i+1}}(EA)(\mathbf{BW}^{(i)})^T(\mathbf{Bd}^{(i)})dx $$

$$= \int_{x_i}^{x_{i+1}}(EA)(\mathbf{W}^{(i)T}\mathbf{B}^T)(\mathbf{Bd}^{(i)})dx= \int_{x_i}^{x_{i+1}}(EA)(\mathbf{W}^{(i)} \cdot \mathbf{B}^T)(\mathbf{Bd}^{(i)})dx$$

$$ \beta=\mathbf{W}^{(i)} \cdot (\int_{x_i}^{x_{i+1}}\mathbf{B}^T(EA)\mathbf{B}dx) d^{(i)}$$

$$ k^{(i)}_{2x2}= \int_{x_i}^{x_{i+1}} \underbrace{\mathbf{B}^T(x)}_{2x1} \underbrace{(EA)(x)}_{1x1} \underbrace{\mathbf{B}  (x)}_{1x2} dx$$

$$B(x)= \begin{bmatrix} 1^{st} Term && \frac{1}{L^i} \end{bmatrix}$$ where $$L^{(i)}= x_{i+1} -x_i $$ (which is the length of element i).

Considering EA is a constant
So the stiffness matrix becomes

$$ \mathbf{k}^{(i)} =\frac{EA}{L^{(i)}} \begin{bmatrix}1&&-1\\-1&&1 \end{bmatrix}$$

Transformation of variable coordinates from x to $$ \tilde{x}$$
$$\tilde{x}:= x-x_i$$ and $$d\tilde{x}=dx$$

$$ k^{(i)}= \int_{\tilde{x}=0}^{\tilde{x}=L^{(i)}} \mathbf{B}^T(\tilde{x}) (EA)(\tilde{x}) \mathbf{B}(\tilde{x}) d\tilde{x}$$

The tapered bar can be drawn and solved using the tilde axis

$$ A(\tilde{x}) = N_1^{(i)}(\tilde{x})A_1 + N_2^{(i)}(\tilde{x})A_2$$

$$ E(\tilde{x}) = N_1^{(i)}(\tilde{x})E_1 + N_2^{(i)}(\tilde{x})E_2$$

$$k^{(i)}_{2x2}$$ derived
$$N_1^{(i)}(\tilde{x}) = N_1^{(i)}(x)$$ which was solve eariler in this report.

$$N_2^{(i)}(\tilde{x})= \frac{\tilde{x}}{L^{(i)}}= \bigg\{ \begin{matrix} 0 && at, \tilde{x}=0 \\ 1 && at, \tilde{x}=L^{(i)} \end{matrix}$$

k(1)
Another way to solve the tapered bar is using the average area and modulus of elasticity.

$$\underbrace{\frac{E(A_1+A_2)}{L}}_{Average Area} \begin{bmatrix}1&&-1 \\1&& 1 \end{bmatrix}= \mathbf{k}^{(i)}$$

Remark
Recall the mean value theorem (MVT) and its relation to the centroid:

MVT : $$\int_{x=a}^{x=b} f(x) dx= f(\bar{x})[b-a]$$

for $$\bar{x}\in[a,b]$$

$$\in$$ symbol means in so $$ a \le \bar{x} \ge b$$

centroid$$= \int_A xdA= \bar{x}\int_A = \bar{x}A$$

$$\int_{a}^{b}f(x)g(x)dx= f(\bar{x})g(\bar{x})[b-a]$$

where $$a \le \bar{x} \ge b$$

but $$f(\bar{x}) \ne \underbrace{\frac{1}{b-a} \int_{a}^{b}f(x)dx}_{Average value of f}$$

and $$g(\bar{x}) \ne \underbrace{\frac{1}{b-a} \int_{a}^{b}g(x)dx}_{Average value of g}$$

Modification of 2-bar truss code to accommodate general k(i)

The following matlab code was used to solve the electric pylon. The functions: PlaneTrussElement, PlaneTrussResults,and NodalSoln were obtained from the Fundamental Finite Element Analysis and Applications book website. Here

The element number for the element with the highest tensile stress is represented by Tbn and equals 81. The highest tensile stress is represented by Ts and equals 9.0511e+005

The element number for the element with the highest compressive stress is represented by Cbn and equals 55. The highest compressive stress is represented by Cs and equals -8.6957e+005

This problem is statically indeterminate. One reaction was given but there are too many reactions remaining to solve using the static equations.

title: User:Eml4500.f08.gravy.mmm/HW6


text:
A course in troll sockpuppets/Fergus MacTroll

Fergus MacTroll, Scottish nationalist.


text: User talk:Uf.team.wiki/HW6

title: User talk:Uf.team.wiki/HW6


text: User talk:Uf.team.wiki/HW7

title: User talk:Uf.team.wiki/HW7


text:
A course in troll sockpuppets/Conspiracy theorist troll

Conspiracy theorist troll is AWESOME. HE KNOWS EVERYTHING ABOUT SECRETS, AN WEARS A FOILHAT TO DETER ALIEN AND GOVERNMENTAL WAVES WHICH CORRUPT YOUR BRAINS. HE ALSO LOVES CAPSLOCK!


text:
User:Eml4500.f08.delta 6.castrillon/3-bar space truss

3-Dimensional Space Truss (MATLAB)
A three-bar space truss is discussed in the course textbook by Asghar Bhatti on page 230. This truss system contains elements with these properties:

a1 = 200; a2 = 600; e = 200000; P = 20000;
 * E$$^{(1)}$$=200,000
 * E$$^{(2)}$$=200,000
 * E$$^{(3)}$$=200,000

 * A$$^{(1)}$$=200 mm
 * A$$^{(2)}$$=200 mm
 * A$$^{(3)}$$=600 mm

The truss system has four nodes, three coplanar and supported by a wall, and belong the the three elements, which meet at the free-moving fourth node. Also, a 20-kN force acts upon the fourth global node in the negative y-direction.

The following MATLAB code shows the colaculation of the system reactions and nodal dofs. It is easy to note that only the fourth global node will showacse any displacement due to the force P.

First defining some functions:

These functions serve to set up each element's stiffness matrix and display the results. The bulk of the code is below, which constructs the stiffness matrix, determines the dof values, and the reactions. The results are also displayed below.

The deformation of this space truss due to the force P and its code is shown below:

This is the truss output from the above code.

x-axis perspective:

y-axis perspective:

z-axis perspective:

(-2,-2,3) perspective:

(1.5,2,2) perspective:

Statics Approach
The three-bar space truss explored above proves to be statically determinate because its reactions were readily able to be determined. Determining these results was possible from the fact that there were enough constraints on the truss, i.e. the fixed nodes on the wall. Since this truss system is statically determinate, then a statics approach can be developed.

This truss system is loaded by -P$$\hat{j}$$.

Defining the vectors from global node 4 to each of the other three global nodes:

 * 41=0.96$$\hat{i}$$+1.92$$\hat{j}$$-2$$\hat{k}$$

 * 42=-.144$$\hat{i}$$+1.44$$\hat{j}$$-2$$\hat{k}$$

 * 43=-2$$\hat{k}$$

And their magnitudes:

 * $$|41|=\sqrt{0.96^2+1.92^2+(-2)^2}=2.9339$$

 * $$|42|=\sqrt{(-1.44)^2+1.44^2+(-2)^2}=2.8543$$

 * $$|43|=\sqrt{(-2)^2}=2$$

Defining the forces along each element (negative for compression) in vector notation:

 * T$$_{41}$$ = $$T_{41}$$ λ$$_{41}$$ =$$T_{41}\frac{\mathbf{T_{41}}}{|41|}$$=$$(0.3272\mathbf{\hat{i}}+0.6544\mathbf{\hat{j}}-0.6817\mathbf{\hat{k}})T_{41}$$

 * T$$_{42}$$ = $$T_{42}$$ λ$$_{42}$$ =$$T_{42}\frac{\mathbf{T_{42}}}{|42|}$$=$$(-0.5045\mathbf{\hat{i}}+0.5045\mathbf{\hat{j}}-0.7007\mathbf{\hat{k}})T_{42}$$

 * T$$_{43}$$ = $$T_{43}$$ λ$$_{43}$$ =$$T_{43}\frac{\mathbf{T_{43}}}{|43|}$$=$$(-\mathbf{\hat{k}})T_{43}$$

Applying the knowledge that the forces will all balance out in all directions, we can solve for T$$_{41}$$, T$$_{42}$$, and T$$_{43}$$:

ΣF$$_x$$:

 * 0.3272$$T_{41}$$-0.5045$$T_{42}$$=0

ΣF$$_y$$:

 * 0.6544$$T_{42}$$+0.5045$$T_{42}$$-P=0

ΣF$$_z$$:

 * -$$T_{43}$$-0.6817$$T_{41}$$-0.7007$$T_{42}$$=0

Three unknowns and three equations allows us to find the reactions along the axis of the elements:

 * $$T_{41}$$=-2,314.533 N (compression)

 * $$T_{42}$$=20,374.691 N (tension)

 * $$T_{43}$$=13,214.269 N (tension)

Notice that these exact reaction values are equal to the values attained with the FEA method above.


text:
User:Eml4500.f08.team.Bengtson/teamHW6

FEM via PVW (Continued)
Prove: $$ u(x_{i+1})=d_{i+1} $$ Citing the previous days notes $$\displaystyle u(x_{i+1})=N_i(x_{i+1})d_i+N_{i+1}(x_{i+1})d_{i+1}=(0)d_i+(1)d_{i+1}=d_{i+1}$$ $$\displaystyle u(x_{i+1})=d_{i+1}$$ The two expressions are equal

Apply same interpolation for w(x), ie

$$ w(x)=N_{i}(x)w_{i}+ N_{i+1}(x)w_{i+1} $$

Element Stiffness matrix for element:

$$ (B)=\int_{x_{i}}^{x_{i+1}}[N_{i}(x)w_{i}+ N_{i+1}(x)w_{i+1}][N_{i}(x)d_{i}+ N_{i+1}(x)d_{i+1}]dx $$ This integral represents the elemental stiffness matrix for the element.

$$ N_{i}=\frac{\partial N_{i}(x) }{\partial x} $$

Likewise for $$ N_{i+1}$$

Note: $$ u(x)=[N_{i}(x)N_{i+1}(x)]\begin{bmatrix} d_{i}\\ d_{i+1} \end{bmatrix} $$

$$N_{i}= [N_{i}(x)N_{i+1}(x)]$$ and is a 1x2 matrix

$$ \frac{\partial N_{i}(x) }{\partial x}=[N'_{i}(x)N'_{i+1}(x)]\begin{bmatrix} d_{i}\\ d_{i+1} \end{bmatrix} $$

$$ B(x)= [N'_{i}(x)N'_{i+1}(x)]$$ and is a 1x2 matrix.

Similarly $$w(x)= N(x)\begin{Bmatrix} w_i\\ w_{i+1} \end{Bmatrix} $$

$$ \frac{\partial \mathbf{w}(x)}{\partial x}= \mathbf{B}(x) \begin{Bmatrix} w_i\\ w_{i+1} \end{Bmatrix} $$

Recall the element defintions These are the axial displacements of the element.

$$ \begin{Bmatrix} w_i\\ w_{i+1} \end{Bmatrix}=\begin{Bmatrix} w_1^i\\ w_2^i \end{Bmatrix}=\mathbf{w}^{(i)} $$

$$ (\boldsymbol{B})=\int_{x_i}^{x_{i+1}}{(\mathbf{B}\mathbf{W}^{(i)})(EA)\mathbf{B}\mathbf{d}^{(i)}}dx $$

$$ =w^{(i)}\bullet (\mathbf{k}^{(i)}\mathbf{d}^{(i)}) $$

From this we want to identify the stiffness matrix.

$$ (\boldsymbol{B})=\int_{x_i}^{x_{i+1}}{(\mathbf{B}\mathbf{W}^{(i)})(EA)(\mathbf{B}\mathbf{W}^{(i)}) \bullet \mathbf{B}\mathbf{d}^{(i)}}dx $$

$$ (\mathbf{B}\mathbf{W}^{(i)}) \bullet \mathbf{B}\mathbf{d}^{(i)}= (\mathbf{B}\mathbf{W}^{(i)})^T \mathbf{B}\mathbf{d}^{(i)}$$

As proven before, the transpose of the first matrix multiplied by the second is equivalent to the dot product of both together. $$(\mathbf{B} \mathbf{W}^{(i)})^T=$$ $$\mathbf{W}^{(i)T} $$ $$\mathbf{B}^T= \mathbf{W}^{(i)}\bullet \mathbf{B}^T$$

$$ B= w^{i}\bullet (\int B^T(x)(EA)\mathbf{B}dx)\mathbf{d}^{(i)} $$

$$ (\boldsymbol{K})^i=\int_{x_i}^{x_{i+1}}{(\mathbf{B}(x))(EA)\mathbf{B}(x)}dx $$

HW $$B(x)= \begin{bmatrix}Hw & \frac{1}{L^{i}}\end{bmatrix} $$

$$ L^{(i)}=x_{i+1}-x_i$$ length of element i

Transference of Var(coordinates) from x to $$\tilde{x}$$

$$ \tilde{x}=x_i-x_o $$

$$ d\tilde{x}=dx $$

title: User:Eml4500.f08.team.Bengtson/teamHW6


title: File:FEABBQNov122008.jpg