User:Eas4200c.f08.aeris.krammer/HW6

 See my comments below. Please don't remove these comment boxes; just add your comment to these comment boxes in case you fix the problems.

Matlab code section updated Eas4200c.f08.aeris.holladay 02:58, 9 December 2008 (UTC)

Copyright Violation/Plagiarism updated Contributing Members section updated Eas4200c.f08.aeris.caceres 06:12, 9 December 2008 (UTC)

Proof that $$ I_{yz} = 0 $$ expanded, Figure 2 added. Eas4200c.f08.aeris.mcrae 06:18, 9 December 2008 (UTC)

Dimensional Analysis
Using dimensional analysis, we are able to provide a strategy for choosing relevant data and forming a relationship between them. Let us take the example of normal stress. The standard equation for normal stress is  σ = F/A  where F is force and A is area. Breaking this equation down into components gives us the following:  [f] = F/L  and [A] = L2 where F is the force, L is the length, and [ ] is "The Dimension Of..." When organized together, we get  [f/A] = (F/L)/L2   =  F/L3 also known as force per volume or a body force.

Another example can be used for strain  (Є) . Strain is calculated as  Є = ΔL/L  where ΔL is the total deformation of the original length and L is the original length. This analysis shows that strain is a unitless parameter. This can be verified through dimensional analysis:  Є = du/dx → ΔL/L  therefore,  [Є] = [du]/[dx] = L/L → 1  which confirms that strain is indeed unitless.

This idea of dimensional analysis has many practical uses. Take for instance Poisson's Ratio  (υ) 



 A = Єxx  B = Єyy   C = σy   D = σx Poisson's Ratio is the ratio of transverse contraction strain to longitudinal extension strain in the direction of stretching force. [Rod Lakes, University of Wisconsin]

Taking the y-component of normal strain at point y and relating it to the x-component of normal strain at point x gives us:  υ = -Єyy/Єxx  therefore,  [υ] = [Єyy]/[Єxx] → 1  Poisson's Ratio is a unitless parameter and this analysis does validate that reasoning.

Torsional Analysis Continued
We can recall that $$[\frac{\partial \sigma _{ij}}{\partial x_{i}}]=\frac{F}{l^{3}}$$, which is a ratio of the force and volume.

As we have determined earlier, the torque is the moment couple that acts about the longitudinal axis of a shaft. Working on the D part of the Roadmap, we can recall the Prandtl stress function $$\displaystyle \phi$$, which plays the role of a potential function. $$\sigma _{yx}=\frac{\partial \phi }{\partial z}$$ $$\sigma _{zx}= -\frac{\partial \phi }{\partial y}$$ $$\displaystyle \sigma _{yx}$$ and $$\displaystyle \sigma _{zx}$$ are components of the gradient of the Prandtl function with respect to $$\displaystyle y$$ and $$\displaystyle z$$. We can also determine that there is a scalar function $$\displaystyle f(x, z, y)$$ and the gradient of this  function is a vector and it can be expressed by  the following equation. $$\triangledown 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} $$ We know that  $$\displaystyle \sigma _{yx}$$ and  $$\displaystyle \sigma _{zx}$$   have nonzero values. We can assume uniform stress through the length of the bar, and having static equilibrium, we can write <p style="text-align:center;">$$\frac{\partial \sigma _{xz} }{\partial z}+\frac{\partial \sigma _{yx} }{\partial y} =0$$

Now, we can now substitute the previously determined  $$\displaystyle \sigma _{yx}$$ and  $$\displaystyle \sigma _{zx}$$   in terms of  $$\displaystyle \phi$$  into this equation: <p style="text-align:center;">$$\frac{\partial }{\partial y}(\frac{\partial \phi }{\partial z}) + \frac{\partial }{\partial z} (-\frac{\partial \phi }{\partial y}) =0$$ This equation can be rearranged and expressed in the following form: <p style="text-align:center;">$$\frac{\partial ^{2}\phi }{\partial y \partial x}=\frac{\partial ^{2}\phi }{\partial z \partial y}$$ Since the function  $$\displaystyle \phi$$  is continuous and smooth, the second derivative is interchangeable. As an example, let's look at the equation below: <p style="text-align:center;">$$\frac{\partial ^{2} \phi}{\partial y \partial z}= \frac{\partial ^{2} \phi}{\partial z \partial y}$$ With the help of the earlier determined strain-displacement relations and stress - strain relations we can write the following relation in terms of the Prandtl function: <p style="text-align:center;">$$\frac{\partial ^{2} \phi}{\partial y ^{2}} + \frac{\partial ^{2} \phi}{\partial z ^{2}} = -2G \theta$$ The left side of this equation is the Laplacian $$\triangledown  ^2 \phi $$

Now, we can derive the stress vector or traction <p style="text-align:center;"> $$\displaystyle t_{3x1}=\sigma _{3x3}n_{3x1}$$ Where $$\displaystyle t$$ is the traction force component, $$\sigma$$ is the stress tensor component, and $$\displaystyle n$$ is the component normal to the surface. Examining the general 2D case first as shown in the figure, let's take a look at the uniform bar with an arbitrarily shaped cross section. Or we can also take a look at a more specific example as shown in the second image.

The stress vector can be determined  on the surface with the help of the normal vector  $$\bar{n}$$. Since we are working in the 2D zy-plane, we can see that $$\displaystyle n_{x}=0$$. Therefore, we can write the following relations <p style="text-align:left;">$$t_{x}=\sigma _{yx}n_{y}+\sigma _{zx}n_{z}=\frac{\partial \phi}{\partial y}n_{y}-\frac{\partial \phi}{\partial z}n_{z}$$ <p style="text-align:left;">$$\displaystyle t_{y}=0 $$ <p style="text-align:left;">$$\displaystyle t_{z}=0$$

Also, $$\displaystyle n_{y}$$ and $$\displaystyle n_{z}$$ can be expressed with the help of the angle the horizontal and vertical axis makes the tangent to the contour line at the examined point. Let us call this angle $$\displaystyle \eta$$. Then <p style="text-align:left;">$$n_{y}= cos \eta =\frac{\partial y}{\partial s}$$ and <p style="text-align:left;">$$n_{z}= sin \eta =-\frac{\partial z}{\partial s}$$ Therefore , <p style="text-align:center;">$$t_{x}=\frac{\partial \phi}{\partial z}\frac{dz}{ds}+ \frac{\partial \phi}{\partial y}\frac{dy}{ds}=\frac{d \phi}{ds}$$

The traction boundary condition is $$\displaystyle t_{x}=0$$ and this can be rewritten as $$t_{x}=0=\frac{d \phi}{ds}$$ From this equation we can conclude that $$\displaystyle \phi$$ is a constant on the surface. Thus, the boundary condition can be expressed as  $$\displaystyle \phi =0$$

Also, working on the cross-section shear stresses, considering a differential area, we can determine the total resultant torque <p style="text-align:left;">$$\displaystyle 2\int \int_{A}^{}{} \phi dxdy$$ There is a relationship between the applied torque T and the rate of twist $$\displaystyle \theta$$ that results from this torque, namely <p style="text-align:left;">$$\displaystyle J=\frac{T}{G \theta}$$ In this equation $$\displaystyle J$$ is called the torsion constant.

Deriving Traction Force Equation cont...
Recall that the equation we are trying to derive is the following. <p style="text-align:center;">$$\displaystyle [t]_{3x1}=[\sigma]_{3x3}[n]_{3x1}$$ Traction force is the stress applied to a lateral surface. For example, stresses due to aerodynamic forces on a wing such as drag or even friction can be a traction force. The goal is to essentially derive the three components of the traction force matrix from the stress tensor and the surface's normal vector. This expression will be obtained by geometrical manipulation. Let's consider the following figure to illustrate an infinitesimally small surface element. <p style="text-align:center;"> In the diagram, ds is the length of the lateral surface, n is the unit normal vector of the lateral surface, and the base and height of this triangular element is dy and dz respectively. Now, in order to relate the traction force to the stress tensor components, the angle of the lateral surface with respect to the coordinate system must be utilized. First, the base and height of the triangle will be written in terms of the length of the lateral surface. Thus, it is simple to see the following by simple geometry. <p style="text-align:center;"> $$\displaystyle dz=dscos(\theta)$$

$$\displaystyle dy=dssin(\theta)$$ Similarly, the normal vector can be written in terms of the angle $$\theta$$. Realize that the magnitude of the normal vector is one, the following can be written. <p style="text-align:center;"> $$||\vec{n}||=1$$ $$\displaystyle n_y=cos(\theta)$$ $$\displaystyle n_z=sin(\theta)$$ To continue with the derivation, the next logical step is to write the equilibrium equations of the forces in the z and y directions. To simplify things, the depth in and out of the page will be assumed to be one (the area of the base and the height will be the base and the height respectively). Although this is the derivation for the 2D case, the process here can be extrapolated with another axis and equations for the 3D case can also be derived. Let's equate forces in the z and y directions. <p style="text-align:center;"> $$\sum {F}_y=0=-\sigma_{yy}dz-\sigma_{yz}dy +t_yd_s$$ Substituting the base and height will the lateral surface length and its respect direction cosine, we get the following equation. $$\displaystyle =-\sigma_{yy}ds n_y-\sigma_{yz}ds n_z +t_yd_s$$ It is possible to divide out the ds terms resulting in... $$\displaystyle t_y= \sigma_{yy}n_y+\sigma_{yz}n_z ... (1)$$ The sum of the forces in the z direction can be found in the below collapsible box. Once the equations 1 and 2 are obtained, this can be arranged into a matrix notation. In book the subscripts are different in that it uses the (x,y) coordinate system. Using cyclic permutation, we have the following equation. <p style="text-align:center;"> $$\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) $$

The above equation is true because $$\sigma_{yz}=\sigma_{zy}$$. Equation 3 is written this way so that the first subscript indicates the row of the equation. In indicial notation, the equation is rewritten as the following. <p style="text-align:center;"> $$\begin{bmatrix} t_1\\ t_2\\ \end{bmatrix} = \begin{bmatrix} \sigma_{11} & \sigma_{12}\\ \sigma_{21} & \sigma_{22}\\ \end{bmatrix} \begin{bmatrix} n_1\\ n_2\\ \end{bmatrix}...(3') $$

This can be projected into the third dimension. Generalizing this equation into the 3D case in indicial notation yields the following. <p style="text-align:center;">$$ \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) $$ The advantage of using the indicial notation is that it can be written in a compact way that still provides all the information of the matrix notation. <p style="text-align:center;"> $$t_i=\sum_{j=1}^3{\sigma_{ij}n_j},i=1,2,3$$ Now, the original equation is derived. <p style="text-align:center;">$$\displaystyle [t_i]_{3x1}=[\sigma_{ij}]_{3x3}[n_i]_{3x1} $$     From the road map, we have the following three equations. <p style="text-align:center;">$$\displaystyle T=2\int\int_{A}^{}{}\phi dA$$ $$\displaystyle T = GJ\theta $$ $$\displaystyle J=\frac{-4}{\nabla^{2}\phi }\int \int _{A}\phi dA$$ Combining these three equations and solving for the Lagrangian of phi, we have the following equation. <p style="text-align:center;">$$ \frac{\partial^2{\phi}}{\partial{y}^2}+\frac{\partial^2{\phi}}{\partial{z}^2}=-2G\theta $$ where $$\phi$$ is constant on the lateral surface of the bar.

Bending and Torsional Analysis Continued
To continue torsional analysis, let us consider two cases of different cross sections. Case 1 is a thin-walled cross section that is closed. An example of this kind of cross section would be a NACA airfoil. It is two lateral surfaces denoted by S with corresponding subscripts as shown in the figure below. Case 2 is a solid cross section with only one lateral surface. Both cross sections are uniform and non circular. <p style="text-align:center;"> In case 2, solid sections with a single contour boundary have phi as an arbitrary constant and thus can be chosen to be zero. This is due to a traction free boundary condition meaning that the traction force t_z is equal to zero. In case 1, using the Prandtl stress function phi, the stress-free boundary condition says that the derivative of phi with respect to the contour length is zero on lateral surfaces o and i. This means that phi on S_o and S_i are different constants. Differing from case 2, these two constants cannot be set to equal zero. Now consider a uniform bar with circular, solid cross section. The Prandtl stress function of this geometry is as follows. <p style="text-align:center;"> <p style="text-align:center;">$$ \phi(y,z)=c\left(\frac{y^2}{a^2}+\frac{z^2}{b^2}-1\right) $$ where a=b=radius and thus constituting a circle. If a did not equal b, the shape of the cross section would be an ellipse. Applying equation 3.19 in the book, it is possible to solve for the factor C. <p style="text-align:center;">$$
 * }

\nabla^2{\phi}=\frac{\partial^2{\phi}}{\partial{y}^2}+\frac{\partial^2{\phi}}{\partial{z}^2}=-2G\theta$$ $$\nabla^2{\left[C(\frac{y^2}{a^2}+\frac{z^2}{a^2}-1)\right]}={C}\nabla^2{\left[(\frac{y^2}{a^2}+\frac{z^2}{a^2}-1)\right]}$$ $$=C\nabla{(\frac{2y}{a^2}+\frac{2z}{a^2})}=\frac{2C}{a^2}\nabla{(y+z)}=\frac{4C}{a^2}$$ Applying the compatibility equation yields the following $$C=-\frac{a^2}{2}G\theta$$

Bars with Circular Cross-Sections
Proving that a bar with a circular cross-section does not warp along its longitudinal axis or z-axis, one must first consider the contour shape of the boundary.

In this case, a circle

<p style="text-align:center;"> represented by the equation... $$x^{2}+y^{2}=r^{2}$$

...is used to express the stress function for the circular bar, with 'C' being defined as...

<p style="text-align:center;"> $$C=-\frac{1}{2}r^{2}G \theta$$

...from the compatibility equation of torsion.

The stress function is expressed as...

<p style="text-align:center;"> $$\phi =C\left(\frac{x^{2}}{r^{2}}+\frac{y^{2}}{r^{2}}-1 \right)$$

Using this formula...

<p style="text-align:center;"> $$T=2\int _{A}\int \phi dx dy$$

...one can calculate the torque in terms of 'C' and the circular contour. The torque can then be expressed as...

<p style="text-align:center;"> $$T=\theta \cdot J \cdot G$$

Acknowledging that shear stress is simply...

<p style="text-align:center;"> $$\tau =\frac{Tr}{J}$$

...the torque is then divided by the polar moment of inertia and multiplied by the radius, or 'X' or 'Y'. It can be seen that the shear stress is equivalent to...

<p style="text-align:center;"> $$\tau_{xz} =\sigma _{xz}=-G\theta y=2C\frac{y}{a^{2}}=\frac{\partial \phi }{\partial y}$$

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

It is seen that...

<p style="text-align:center;"> $$\sigma _{xz}=-G \cdot \theta \cdot y$$

$$\sigma _{yz}=G \cdot \theta \cdot x$$

can be substituted into the strain-displacement relations...

<p style="text-align:center;"> $$\gamma _{xz}=\frac{\sigma_{xz}}{G}=\frac{\partial u_{z}}{\partial x}-\theta y$$

$$\gamma _{yz}=\frac{\sigma_{yz}}{G}=\frac{\partial u_{z}}{\partial y}+\theta x$$

...simplifying to...

<p style="text-align:center;"> $$\frac{\partial u_{z}}{\partial x}=0$$

$$\frac{\partial u_{z}}{\partial y}=0$$

...which results in...

<p style="text-align:center;"> $$u_{z}=0$$ .


 * }

<div style="width: 80%; margin-left: auto; margin-right: auto; padding: 4px; border: 2px solid #FF0000; background-color: #FFDDDD; text-align: left;"> Error: See my comments in Team VQCrew (with annotation). Eml4500.f08 14:04, 23 November 2008 (UTC)

Shear Flow in Open, Thin-Walled Sections


General equation for asymmetric cross sections:
 * $$\ \int_{A_{s}}^{}{ \frac{d \sigma_{xx}}{dx}dA} = -q_{s}$$

Symmetric Case:


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


 * $$\ q_{s} = - \frac{V_{z}}{I_{y}} \int_{A_{s}}^{}{zdA} = \frac{-V_{z}Q_{y}}{I_{y}} $$        (eq. 5.2)


 * $$\ Q_{y} = \int_{A_{s}}^{}{zdA} = z_{c}A_{s} $$

Asymmetric Case:


 * $$\ \sigma _{xx} = \left(k_{y}M_{z} - k_{yz}M_{y} \right)y + \left(k_{z}M_{y} - k_{yz}M_{z} \right)z $$         (eq. 5.3)

with


 * $$\ k_{y} = \frac{I_{y}}{D} $$

and


 * $$\ k_{yz} = \frac{I_{yz}}{D} $$

and


 * $$\ k_{z} = \frac{I_{z}}{D} $$

This can be also written in matrix form:


 * $$\ \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} $$


 * $$\ \sigma _{xx} = \begin{bmatrix} y & z \end{bmatrix}\begin{bmatrix} k_{y}M_{z} & -k_{yz}M_{y}\\ -k_{yz}M_{z} & k_{z}M_{y} \end{bmatrix} $$

We can also flop the order of variables, and come up with the preferred form of this equation:


 * $$\ \sigma _{xx} = \begin{bmatrix} z & y \end{bmatrix}\begin{bmatrix} k_{z} & -k_{yz}\\ -k_{yz} & k_{y} \end{bmatrix} \begin{Bmatrix} M_{y} \\ M_{z} \end{Bmatrix} $$

We particularize to symmetric cross section, and find that $$\ I_{yz} = 0 $$. Proof of this is shown in the tab below.

Consider $$\ M_{z} = 0 $$

We know the following relations.


 * $$\ I_{yz} = 0 $$
 * $$\ D = I_{y}I_{z} $$
 * $$\ k_{yz} = 0 $$
 * $$\ k_{z} = \frac{1}{I_{y}} $$

Combining these equations we find that $$\ \sigma_{xx} = \frac{M_{y}z}{I_{y}} $$.

Asymmetric Case


 * $$\ q_{s} = -\left( k_{y}V_{y} - k_{yz}V_{z} \right) Q_{z} - \left( k_{z}V_{z} - k_{yz}V_{y} \right) Q_{y} $$          (eq. 5.5)


 * $$\ Q_{z} = \int_{A_{s}}^{}{y*dA} $$


 * $$\ Q_{y} = \int_{A_{s}}^{}{z*dA} $$

For stringer web sections, if the thickness (t) of the skin and spar webs is very small we can neglect the computations of Iy, Iz, Iyz, Qy, and Qz. We use only the areas of the stringers.

HW6 MatLAB Problem
For this homework assignment we take a closer look at the effects of plate buckling. We begin by modifying the NACA airfoil code to calculate the stresses along the top and bottom skin panels. Using the equation for bending stress we calculate the stresses from point B to F and then from point E to H. The results are diplayed below. From the graph we can see that the top panel is under compression while the bottom panel is under tension. Also, the maximum absolute stress occurs near point B as determined last time.

<div style="width: 80%; margin-left: auto; margin-right: auto; padding: 4px; border: 2px solid #FF0000; background-color: #FFDDDD; text-align: left;"> The figures you plotted were not the perspective view of the buckling mode shapes; see Team Aero. Eml4500.f08 15:42, 23 November 2008 (UTC)

The graph above is a view along the edge of a rectangular panel of dimensions 1 m by 1.5 m. From the graph we can can see the effect of the m in the number of half waves. Another way to determine m is using the period of the function. For the second graph (m=2), one complete cycle is shown so the period is $$a$$ whereas the first graph has period $$2a$$. Thus respectively, $$m=1$$ and $$m=2$$. Plotting the coefficient k_c versus the aspect ratio a/b for different values of m it becomes evident that the number of half waves affects the width of k_c at minimum values. When plotting the critical stress values for our airfoil panel, we see that the minimum stress value occurs when $$a/b=1$$. By looking at the equation for critical stress we can see why this is true. $$\sigma_{cr}=k_{c}*\frac{\pi^2*D}{b^2*h}$$ Since $$b^2$$ is in the denominator it will greatly increase the critical stress at lower values, hence the steep increase to the left of $$a/b=1$$. From the previous graph for $$m=1$$ we see that k_c is minimized at $$a/b=1$$, hence the minumum in critical stress at $$a/b=1$$.

Next we look at the clamped boundary condition for plate buckling.

$$  \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) $$

Looking at the edge (y = b) of the panel we do not see any displacement in the z-direction for the clamped boundary condition. Examining the equation for Ψ we see that at y = b, the second term becomes $$cos(2\pi)$$. Therefore, the displacment along the y = 0 and y = b edge of the plate will be zero. However, if we look at the profile at y/2 we can see the maximum displacement as displayed below.

Contributions
Jesse W. McRae Eas4200c.f08.aeris.mcrae 19:43, 21 November 2008 (UTC) Josh Holladay Eas4200c.f08.aeris.holladay 19:57, 21 November 2008 (UTC) Anett Krammer Eas4200c.f08.aeris.krammer 20:34, 21 November 2008 (UTC) Jin Yu Guan Eas4200c.f08.aeris.guan 21:34, 21 November 2008 (UTC) Nelson Caceres Eas4200c.f08.aeris.caceres 21:45, 21 November 2008 (UTC)