text:
User:Eml5526.s11.team5/sub3

= 3.1 Weighted Residue Form as BC on Polynomial Base Function (k = 1) =

Given

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

$$ \omega = ]0,1[ \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 3.1.1)
 * }

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

$$ a_2 = 2 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 3.1.2)
 * }

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

$$ f = 3 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 3.1.3)
 * }

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

$$ \frac{\partial{u}}{\partial{t}}(x,t) = 0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 3.1.4)
 * }

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

$$ \Gamma_g = 1 $$, $$ g = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle\!
 * }

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

$$ u(x=1) = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 3.1.5)
 * }

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

$$ \Gamma_h = 0 $$, $$ h = 4 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ -\frac{du^h}{dx}(0) = 4 \quad $$ $$
 *  $$ \displaystyle (Eq. 3.1.6)
 * }

Find
 Find approximate solution $$ u^n \quad $$ and compare it to exact one

Consider base function $$ b_j(x), j = 0, 1, 2, 3,... n \quad $$

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

$$ b_j(x) = (x+k)^j \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.7)
 * }

,where $$ k \quad $$ is an arbitrary constant, which in this exercise it chosen to be $$ 1 \quad $$

(1) Let $$ n = 2 \Rightarrow ndof = n+1 = 3 \quad $$

Generate variables considered for $$ n = 2 \quad $$ case

(2) Find 2 eqs that enforce b.c.'s for $$ u^h(x) = \sum_{j -0}^{n} d_j b_j(x) $$

(3) Find 1 more eq. to solve for $$ \underline{d} = d_j (j = 0, 1, 2) $$ by projecting the residue $$ \underline{p}(u^h) $$ on a basis function $$ b_k(x) \quad $$ with $$ k = 0, 1, 2 \quad $$. The additional equation must be lineally independent from the above 2 equations in part (2).

(4) Display 3 equations in matrix form

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

$$ \underline{K}\underline{d} = \underline{F} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.8)
 * }

Observe if $$ \underline{K} $$ is symmetric

(5) Solve for $$ \underline{d} \quad $$

(6) Construct $$ u_n^k(x) $$ and plot $$ u_n^k(x) $$ vs. $$ ux) $$

(7) Repeat steps (1) to (6) for

 * (7.1) $$ n = 4 \quad $$
 * (7.2) $$ n = 6 \quad $$

(8) Compute $$ u_n^h (x = 0.5) \quad $$ for $$ n = 2, 4, 6 \quad $$

error

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

$$ e_n(0.5) = u_n(0.5) - u_n^h(0.5) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.9)
 * }

Plot $$ e_n(0.5) \quad $$ vs. $$ n \quad $$

Solution
(1)

For $$ n = 2 \quad $$

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

$$ \underline{b} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} b_0     & b_1 & b_2 \end{bmatrix}

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

From base function expression $$ Eq. 3.3.7 \quad $$

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

$$ \underline{b} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} (x+1)^0     & (x+1)^1 & (x+1)^2 \end{bmatrix}

=

\begin{bmatrix} 1     & x+1 & (x+1)^2 \end{bmatrix}

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

Unknown function

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

$$ \underline{d} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} d_0     & d_1 & d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle (Eq. 3.1.12)
 * }

Approximate solution

(2)

at B.C. $$ u(1) = 0 \quad $$

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

$$ u^h = \sum_{i=0}^{n=2} d_ib_i $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.13)
 * }

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

$$ u^h(x) = d_0(1)   + d_1(x+1) +d_2(x+1)^2   \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.14)
 * }

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

$$ u^h(1) = 0 = d_0(1)   + d_1(1+1) +d_2(1+1)^2 =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1     & 2 & 4      \end{bmatrix}

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle\!
 * }

at B.C $$ \frac{du}{dx}(0) = 4 $$

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

$$ (u')^h(x) = \frac{du}{dx} = \sum_{i=0}^{n=2} d_ib_i' $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.15)
 * }

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

$$ \underline{b}' =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0     & 1 & 2(x+1) \end{bmatrix}

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

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

$$ (u')^h(0) = - 4 = d_00 + d_11 + d_22(0+1) =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0     & 1 & 2     \end{bmatrix}

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle\!
 * }

(3)

By projecting the residue part of a function on the bases, 3rd equation can be generated:

From the definition of elastodynamics case of ODE function $$ P_u(x) \quad $$

Rewriting $$ Eq. 3.1.7 \quad $$

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

$$ b_j(x) = (x+k)^j \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.7)
 * }

In this, more specific, "self adjoin" case

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

$$ P(u^h) :=[a_2u']' + a_0u - g = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.17)
 * }

noting that

$$ g = 0 \Rightarrow \quad $$ no forcing function (steady state)

$$ a_0u = f = 3 \quad $$ (given)

$$ a_2 = 2 \quad $$ (given)

$$ P(u^h) :=[2u']' + 3 = 0 \Rightarrow \quad $$ since $$ a_2 \quad $$ is given as constant 2 (not function of "x"), it can be taken out of bracketing differential, resulting in

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

 * style="width:95%" |
 * style="width:95%" |

$$ P(u^h) :=2u'' + 3 = 0= 2\frac{d^2u}{dx} +3 = 0 $$ $$
 *
 $$ \displaystyle (Eq. 3.1.18)
 * }

From $$ Eq. 3.1.13 \quad $$, $$ u^h \quad $$ was defined as

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

$$ u^h = \sum_{i=0}^{n=2} d_ib_i = u^h(x) = d_0(1)   + d_1(x+1) +d_2(x+1)^2   \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.19)
 * }

Taking $$ 2^{nd} \quad $$ derivative with respect to "x"

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

$$ u^{h''} = d_00 + d_10 + 2d_2 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.20)
 * }

Therefore

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

$$ P(u^h) := 2(d_00 + d_10 + 2d_2)) +3 = 4d_2 + 3\quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.21)
 * }

Afterwords, noting that

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

$$ \int_\omega w P(u^h) dx := 0 \Rightarrow w_i = b_i \Rightarrow \int_\omega b_i P(u^h) dx := 0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.22)
 * }

which is also noted in $$ 10-4, Eq. (2) \quad $$

and substituting for $$ P(u^h) \quad $$ equation above, we get

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

 * style="width:95%" |
 * style="width:95%" |

$$ \int_\omega b_i (4d_2 + 3) dx := 0 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.1.23)
 * }

recalling that $$ \underline{b} =

\begin{bmatrix} 1  \\ (x+1) \\ (x+1)^2 \end{bmatrix} $$

and taking any of the equations as the last needed equation, for example $$ b_1 \quad $$, we get

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

$$ \int_\omega
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1  \\ (x+1) \\ (x+1)^2 \end{bmatrix}

[0d_0 + 0d_1 + 4d_2] dx := \int_\omega (-3

\begin{bmatrix} 1  \\ (x+1) \\ (x+1)^2 \end{bmatrix} )dx \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.1.25)
 * }

Doing symbolic integration and evaluating interval on $$ \omega = ]0,1[ \quad $$ we get

The outputs are as follow:

$$ \int_\omega

\begin{bmatrix} 1  \\ (x+1) \\ (x+1)^2 \end{bmatrix}

[0d_0 + 0d_1 + 4d_2] dx = $$ 0        0    4.0000         0         0    6.0000         0         0    9.3333

$$ = \int_\omega (-3

\begin{bmatrix} 1  \\ (x+1) \\ (x+1)^2 \end{bmatrix} )dx \quad $$

-3.0000  -4.5000   -7.0000

*********************************************

(4)

Assembling the matrix $$ \underline{K} \underline{d} = \underline{F} \quad $$ from 2 B.C. constrains and 1 more equation generated by weighting factor method:

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

$$ \begin{bmatrix} 1     & (x+1) & (x+1)^2 \\ 0 &   1  & 2(x+1) \\ 0 &   0  & 4     \end{bmatrix}
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix}

=

\begin{bmatrix} 0\\ 4 \\ -3   \end{bmatrix} $$ $$
 *
 $$ \displaystyle (Eq. 3.1.26)
 * }

For "K" matrix: 1    2     4     0     1     2     0     0     4

For "F" matrix:

0   -4    -3

(5)

Solving for $$ \underline{d} = \underline{K}^{-1} F\quad $$ we get

d = 8.0000  -2.5000   -0.7500

(6)

Constructed $$ u^h \quad $$ matrix on matlab by using formula $$ u^h = \sum_{i=0}^{n=2} d_ib_i $$

8.0 - 2.5*x - 0.75*(k + x)^2 - 2.5*k

Exact solution:

Given a general expression for 2nd Order ODE

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

$$ a_2(x)y'' + a_1(x)y + a_0(x)y = g(x) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.27)
 * }

and taking a more specific case ("Self adjoin") for this specific problem, we have

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

$$ [a_2(x)y']' + a_0(x)y = g(x) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.17)
 * }

$$ g(x) = 0 \Rightarrow \quad $$ no forcing function (steady state)

$$ a_0(x)y = f(x) = 3 \quad $$ (given)

$$ a_2 = 2 \quad $$ (given)

Therefore $$ Eq. 3.1.17 \quad $$ becomes

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

 * style="width:95%" |
 * style="width:95%" |

$$ [2y']' + 3 = 0 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.1.28)
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ \frac{du}{dx}(\frac{du}{dx}2) + 3 = 0 $$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ \int \frac{du}{dx}2 = \int -3 dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \frac{du}{dx}2 = -3x + c_1' \Rightarrow \frac{du}{dx} = -1.5x + c_1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.29)
 * }

at B.C of $$ -\frac{du}{dx}(0) = 4 \Rightarrow $$

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

 * style="width:95%" |
 * style="width:95%" |

$$ -4 = -1.5(0) + c_1 \Rightarrow c_1 = -4$$ $$
 *
 $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ \int {du} = \int (-1.5x - 4) dx $$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ u = \frac{-1.5}{2}x^2 - 4x + c_2 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.30)
 * }

at B.C of $$ u(1) = 0 \Rightarrow $$

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

 * style="width:95%" |
 * style="width:95%" |

$$ 0 = \frac{-1.5}{2}(1)^2 - 4(1) + c_2 \Rightarrow c_2 = -4.75$$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ u(x) = -0.75x^2 - 4x - 4.75 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.1.31)
 * }

Circle view

Plain line view

Zoomed in view

So taken in account the exact solution formula, different plot were generated for different n values (n = 2, 4, 6)

Error came out to be:

x = 0.5

ex = n=2       n=4          n=6 0    7.0610e-014  -2.6972e-009

Note, that the 2 member bases approximated the given problem more exact then n = 4 or n = 6. To be exact, it approximated exactly. This was unusual, it might have had been a case that the point around which the error was calculated just happen to lay at that point. So we checked for x = 1, x = 0.25 and x = 0.75. Here are the results:

x = 1

ex = n=2       n=4          n=6 0    1.1369e-013  -3.5798e-009

x = 0.25

ex = n=2       n=4          n=6 0    5.8176e-014  -1.9985e-009

x = 0.75

ex = n=2       n=4          n=6 0    8.8818e-014  -3.2879e-009

What results basically advocated was that n =2 exactly approximated 2nd order ODE, while for this base, increasing number of members 'n' produce relatively greater error (see graphs below).

Netherless, from previous experiences such as convergence of solutions that are seen below and re-use of the earlier proven to work code, it is clear that for these bases $$ b_j = (x+k)^j \quad $$ higher order terms induce relatively higher error then lower order terms.

Approximation and error outputs for bases $$ b_j(x) = (x+k)^j \quad $$

Here is an interesting observation:

From the code attached (line 153 and 154), you may see from error analysis that n=2 and exact value should be the same. Netherless, the zoomed in plot gives different results. It basically says that other approximations are more exact then n = 2. On the other hand, I also added an error (0.5) points of of numerical (black dot) and symbolic (red circle) values.

approxzz = plot(0.5,subs(uhConstest(1),x,0.5),'k.', 'MarkerSize', 15); truezz = plot(0.5,subs(u,x,0.5),'ro', 'MarkerSize', 15);

It clearly shows that these 2 points do coincide and they coincide on the approximated (n = 2) line instead of exact curve. This clearly shows that symbolic value is not as exact as it might have had been presumed it was.

Approximation zoomed in view 251 equidistant point were taken for non-symbolic plot of n=2,4,6

Approximation zoomed in view 51 equidistant point were taken for non-symbolic plot of n=2,4,6

The 2 zoomed in figures above shows how different a number of x values (nodes) affect quality of a plot ('exact' solution curve approaching the exact solution point).

--Egm6341.s11.team3.JA 02:32, 10 February 2011 (UTC)

= 3.2 For the spring mass system given in Exercise 2.1 in Fish and Belytschko =

Find
a) Number the elements and nodes.

b) Assemble the Global Stiffness Matrix and Force Matrix.

c) Partition the system and solve for the Nodal Displacements.

d) Compute the Reaction Forces.

$$ \mathbf{Part(d) : Computing \quad of \quad Reaction \quad Forces. } $$
Lokeshdahiya 15:11, 16 February 2011 (UTC)

= 3.3 Finding the nodal displacements, stresses and reaction forces for a 2-D truss =

Given

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

Considering the truss structure shown in Figure 3.3.1; Nodes A and B are fixed. A force equal to 10N acts in the positive x-direction at node C. Coordinates of joints are given in meters. Young's modulus is $$E=10^{11}Pa\quad$$ and the cross-sectional area for all bars are $$A=2*10^{-2}m^{2}\quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

Find

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

A) Number the elements and nodes.
 * style="width:95%" |
 * style="width:95%" |

B) Assemble the global stiffness and force matrix.

C) Partition the system and solve for the nodal displacements.

D) Compute the stresses and reactions. $$
 *
 $$ \displaystyle\!
 * }

Solution

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

A) Numbering the elements and nodes of the truss structure shown in Figure 3.3.1 is an elementary procedure. The arrangement of elements and nodes in a truss structure is arbitrary, therefore their numbering may be accomplished without any calculation, as seen in Figure 3.3.2 below. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Fig 3.3.1)
 * }

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

$$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Fig 3.3.2)
 * }

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

B) To assemble the global and stiffness force matrix, the truss in Figure 3.3.2 may be further broken down into it's separate constituent elements. Figure 3.3.3 below shows each local element and global coordinate system.  $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Fig 3.3.3)
 * }

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

Element 1 has been numbered with global nodes 1 and 4. Shown in Figure 3.3.3 it is positioned at an angle $$\phi^{(1)}=90^{o}\quad$$, with respect to the x-axis. With $$cos(90^{o})=0, \ \ sin(90^{o})=1, \ \ l^{(1)}=1m, \ \ k^{(1)}=\frac{A^{(1)}E^{(1)}}{l^{(1)}}=\frac{AE}{1m}\quad$$, the element 1 local stiffness matrix $$\textbf{K}^{(1)}\quad$$ may be solved: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{K}^{(1)}=\frac{AE}{1m} \begin{matrix}
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0 & 0 & 0 & 0\\ 0 & 1 & 0 & -1\\ 0 & 0 & 0 & 0\\ 0 & -1 & 0 & 1 \end{bmatrix}

\begin{matrix} \\ (1)\\ \\ (4)\\ \\ \end{matrix}

\\

\begin{matrix} &(1) & & (4)  & \end{matrix}

\end{matrix}$$ $$
 *
 $$ \displaystyle\!
 * }

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

The local stiffness matrices for all elements may be derived from the following equation, courtesy of "A First Course in Finite Elements," by Jacob Fish and Ted Belytschko, 2009, pg. 30: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{K}^{(e)}=k^{e}\begin{bmatrix} cos^{2}\phi^{e} &cos\phi^{e}sin\phi^{e} &-cos^{2}\phi^{e}  &-cos\phi^{e}sin\phi^{e} \\ cos\phi^{e}sin\phi^{e} &sin^{2}\phi^{e} &-cos\phi^{e}sin\phi^{e}  &-sin^{2}\phi^{e} \\ -cos^{2}\phi^{e}&-cos\phi^{e}sin\phi^{e} &cos^{2}\phi^{e}  &cos\phi^{e}sin\phi^{e} \\ -cos\phi^{e}sin\phi^{e} &-sin^{2}\phi^{e} &cos\phi^{e}sin\phi^{e}  &sin^{2}\phi^{e} \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

In a similar fashion, element 2 has been numbered with global nodes 4 and 3. Shown in Figure 3.3.3 it is positioned at an angle $$\phi^{(2)}=0^{o}\quad$$, with respect to the x-axis. With $$cos(0^{o})=1, \ \ sin(0^{o})=0, \ \ l^{(2)}=1m, \ \ k^{(2)}=\frac{A^{(2)}E^{(2)}}{l^{(2)}}=\frac{AE}{1m}\quad$$, the element 2 local stiffness matrix $$\textbf{K}^{(2)}\quad$$ may be solved: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{K}^{(2)}=\frac{AE}{1m} \begin{matrix}
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1 & 0 & -1 & 0\\ 0 & 0 & 0 & 0\\ -1 & 0 & 1 & 0\\ 0 & 0 & 0 & 0 \end{bmatrix}

\begin{matrix} \\ (3)\\ \\ (4)\\ \\ \end{matrix}

\\

\begin{matrix} &(3) & & (4)  & \end{matrix}

\end{matrix}$$ $$
 *
 $$ \displaystyle\!
 * }

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

Also, element 3 has been numbered with global nodes 2 and 4. Shown in Figure 3.3.3 it is positioned at an angle $$\phi^{(3)}=135^{o}\quad$$, with respect to the x-axis. With $$cos(135^{o})=\frac{-1}{\sqrt 2}, \ \ sin(135^{o})=\frac{1}{\sqrt 2}, \ \ l^{(3)}=\sqrt 2m, \ \ k^{(3)}=\frac{A^{(3)}E^{(3)}}{l^{(3)}}=\frac{AE}{\sqrt 2m}\quad$$, the element 3 local stiffness matrix $$\textbf{K}^{(3)}\quad$$ may be solved: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{K}^{(3)}=\frac{AE}{\sqrt 2m} \begin{matrix}
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2}\\ -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} & -\frac{1}{2}\\ -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} & -\frac{1}{2}\\ \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2} \end{bmatrix}

\begin{matrix} \\ (2)\\ \\ (4)\\ \\ \end{matrix}

\\

\begin{matrix} &(2) & & (4)  & \end{matrix}

\end{matrix}$$ $$
 *
 $$ \displaystyle\!
 * }

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

Lastly, element 4 has been numbered with global nodes 2 and 3. Shown in Figure 3.3.3 it is positioned at an angle $$\phi^{(4)}=90^{o}\quad$$, with respect to the x-axis. With $$cos(90^{o})=0, \ \ sin(90^{o})=1, \ \ l^{(4)}=1m, \ \ k^{(4)}=\frac{A^{(4)}E^{(4)}}{l^{(4)}}=\frac{AE}{1m}\quad$$, the element 4 local stiffness matrix $$\textbf{K}^{(4)}\quad$$ may be solved: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{K}^{(4)}=\frac{AE}{1m} \begin{matrix}
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0 & 0 & 0 & 0\\ 0 & 1 & 0 & -1\\ 0 & 0 & 0 & 0\\ 0 & -1 & 0 & 1 \end{bmatrix}

\begin{matrix} \\ (2)\\ \\ (3)\\ \\ \end{matrix}

\\

\begin{matrix} &(2) & & (3)  & \end{matrix}

\end{matrix}$$ $$
 *
 $$ \displaystyle\!
 * }

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

Direct Assembly of the local stiffness matrices yields the global stiffness matrix and the corresponding global displacement, force and reaction matrices: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

$$ K=\frac{AE}{1m} \begin{matrix}

\begin{bmatrix} 0 &0 &0 &0 &0 &0 &0 &0\\ 0 &1 &0 &0 &0 &0 &0 &-1 \\ 0 &0 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & 0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \\ 0 &0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2}+1 & 0 & -1 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} \\ 0 & 0 & 0 & 0 & 1 & 0 & -1 & 0\\ 0 & 0 & 0 & -1 & 0 & 1 & 0 & 0\\ 0 &0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & -1 & 0 & 1+\frac{1}{2} & -\frac{1}{2} \\ 0 &-1 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & 0 & -\frac{1}{2} & 1+\frac{1}{2} \\ \end{bmatrix}

\begin{matrix} \\ (1)\\ \\ (2)\\ \\ (3)\\ \\ (4)\\ \\ \end{matrix}

\\

\begin{matrix} (1)& & &(2)& & &(3)& & &(4) \end{matrix}

\end{matrix} $$

$$ d=\begin{bmatrix} 0\\ 0\\ 0\\ 0\\ u_{3x}\\ u_{3y}\\ u_{4x}\\ u_{4y} \end{bmatrix},

f=\begin{bmatrix} 0\\ 0\\ 0\\ 0\\ 10\\ 0\\ 0\\ 0 \end{bmatrix},

r=\begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y}\\ 0\\ 0\\ 0\\ 0 \end{bmatrix} $$

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

C) The global system may now be partitioned after four rows and four columns, and the nodal displacements may then be solved: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$
 *
 $$ \displaystyle\!
 * }

$$ K_{F}= \begin{bmatrix} 1& 0 & -1 & 0\\ 0 & 1 & 0&0\\ -1 & 0 & \frac{3}{2} & -\frac{1}{2}\\ 0 & 0 & -\frac{1}{2} & \frac{3}{2}\\ \end{bmatrix},

K_{EF}= \begin{bmatrix} 0& 0 & 0 & 0\\ 0& 0 & 0 & -1\\ 0 & 0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2}\\ 0 & -1 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} \end{bmatrix},

d_{F}= \begin{bmatrix} u_{3x}\\ u_{3y}\\ u_{4x}\\ u_{4y} \end{bmatrix},

\bar d_{E}= \begin{bmatrix} \bar u_{1x}=0\\ \bar u_{1y}=0\\ \bar u_{2x}=0\\ \bar u_{2y}=0 \end{bmatrix},

f_{F}= \begin{bmatrix} 10\\ 0\\ 0\\ 0 \end{bmatrix},

r_{E}= \begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y} \end{bmatrix}, $$

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

The unknown displacement matrix is found from the solution of the reduced system of equations below: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\frac{AE}{1m}\begin{bmatrix} 1 &0 &-1  &0 \\ 0 &1  &0  &0 \\ -1 &0  &\frac{3}{2}  &-\frac{1}{2} \\ 0&0 &-\frac{1}{2}  & \frac{3}{2} \end{bmatrix}\begin{bmatrix} u_{3x}\\ u_{3y}\\ u_{4x}\\ u_{4y} \end{bmatrix}=\begin{bmatrix} 10\\ 0\\ 0\\ 0 \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$
 * style="width:95%" |
 * style="width:95%" |

\therefore \begin{bmatrix} u_{3x}\\ u_{3y}\\ u_{4x}\\ u_{4y} \end{bmatrix}=\frac{1m}{AE}\begin{bmatrix} 4 &0 &3  &1 \\ 0 &1  &0  &0 \\ 3 &0  &3  &1 \\ 1&0  &1  & 1 \end{bmatrix} \begin{bmatrix} 10\\ 0\\ 0\\ 0 \end{bmatrix}=\frac{1m}{AE}\begin{bmatrix} 40\\ 0\\ 30\\ 10 \end{bmatrix}=\begin{bmatrix} 2.0\\ 0.0\\ 1.5\\ 0.5 \end{bmatrix}1*10^{-8}m$$ $$
 *
 $$ \displaystyle\!
 * }

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

D) Now in computing the reactions and the stresses: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$r_{E}= \begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y} \end{bmatrix}=K_E\bar d_E+K_{EF}d_F=\begin{bmatrix} 0 &0 &0  &0 \\ 0 &0  &0  &-1 \\ 0 &0  &-\frac{1}{2\sqrt2}  &\frac{1}{2\sqrt2} \\ 0 &-1 &\frac{1}{2\sqrt2}  &-\frac{1}{2\sqrt2} \end{bmatrix}\cdot \frac{1m}{AE}\begin{bmatrix} 40\\ 0\\ 30\\ 10 \end{bmatrix}=\frac{1m}{AE}\begin{bmatrix} 0\\ -10\\ \frac{-10}{\sqrt2}\\ \frac{10}{\sqrt2} \end{bmatrix}=\begin{bmatrix} 0\\ -5\\ -3.53553\\ 3.535534 \end{bmatrix}1*10^{-9}N$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Lastly, in computing each truss constituent element stress, the following equation is implemented (ref: Fish and Belytschko P30, Eq (2.47)): $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \sigma^e = E^e \frac{u_{2x}^{'e}-u_{1x}^{'e}}{l^{e}} = \frac{E^e}{l^e} \begin{bmatrix} -1 & 0 & 1 & 0 \end{bmatrix} \begin{bmatrix} u_{1x}^{'e}\\ u_{1y}^{'e}\\ u_{2x}^{'e}\\ u_{2y}^{'e} \end{bmatrix} =\frac{E^e}{l^e} \begin{bmatrix} -\cos \phi^e & -\sin \phi^e & \cos \phi^e & \sin \phi^e \end{bmatrix} \mathbf d^e $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

For element 1, $$\phi^{(1)}=90^{o}\quad$$, with respect to the x-axis. ($$cos(90^{o})=0, \ \ sin(90^{o})=1$$), and $$\textbf{d}^{(1)}\quad$$ may be equated to: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{d}^{1}=\begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{4x}\\ u_{4y} \end{bmatrix}=\frac{1m}{AE}\begin{bmatrix} 0\\ 0\\ 30\\ 10 \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(1)}=[0\ -1\  0\  1]\begin{bmatrix} 0\\ 0\\ 30\\ 10 \end{bmatrix}\frac{1}{A}=\frac{10}{A}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(1)}=500.00 \ N$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

For element 2, $$\phi^{(2)}=0^{o}\quad$$, with respect to the x-axis. ($$cos(0^{o})=1, \ \ sin(0^{o})=0$$), and $$\textbf{d}^{(2)}\quad$$ may be equated to: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{d}^{2}=\begin{bmatrix} u_{3x}\\ u_{3y}\\ u_{4x}\\ u_{4y} \end{bmatrix}=\frac{1m}{AE}\begin{bmatrix} 40\\ 0\\ 30\\ 10 \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(2)}=[-1\ 0\  1\  0]\begin{bmatrix} 40\\ 0\\ 30\\ 10 \end{bmatrix}\frac{1}{A}=\frac{-10}{A}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(2)}=500.00 \ N, (Positive value due to orientation of element wrt F)$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

For element 3, $$\phi^{(3)}=135^{o}\quad$$, with respect to the x-axis. ($$cos(135^{o})=-\frac{1}{\sqrt 2}, \ \ sin(135^{o})=\frac{1}{\sqrt 2}$$), and $$\textbf{d}^{(3)}\quad$$ may be equated to: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{d}^{3}=\begin{bmatrix} u_{2x}\\ u_{2y}\\ u_{4x}\\ u_{4y} \end{bmatrix}=\frac{1m}{AE}\begin{bmatrix} 0\\ 0\\ 30\\ 10 \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(3)}=\frac{1}{\sqrt2}[\frac{1}{\sqrt2}\ -\frac{1}{\sqrt2}\  -\frac{1}{\sqrt2}\  \frac{1}{\sqrt2}]\begin{bmatrix} 0\\ 0\\ 30\\ 10 \end{bmatrix}\frac{1}{A}=\frac{-20}{A}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(3)}=-1000.00 \ N$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

For element 4, $$\phi^{(4)}=90^{o}\quad$$, with respect to the x-axis. ($$cos(90^{o})=0, \ \ sin(90^{o})=1$$), and $$\textbf{d}^{(4)}\quad$$ may be equated to: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\textbf{d}^{4}=\begin{bmatrix} u_{2x}\\ u_{2y}\\ u_{3x}\\ u_{3y} \end{bmatrix}=\frac{1m}{AE}\begin{bmatrix} 0\\ 0\\ 40\\ 0 \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(4)}=[0\ -1\  0\  1]\begin{bmatrix} 0\\ 0\\ 40\\ 0 \end{bmatrix}\frac{1}{A}=0$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore \sigma^{(4)}=0.00 \ N$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

= 3.4 Non-symmetry of Weighted Residue Form =

Given

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

$$ b_{i(x)}=cos ix \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq.3.4.1)
 * }

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

$$ b_j(x)=cos jx \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.2)
 * }

 and

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

 i≠j
 $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

Find

 Show α ≠ β

Solution

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

$$\alpha = b_i[a^'_{2}b^'_j + a_{2}b''_j]$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.3)
 * }

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

$$\beta = b_j[a^'_{2}b^'_i + a_{2}b''_i]$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.4)
 * }

 Take the first and second derivative of Eq.3.4.1

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

$$ b'_i=-isin(ix) \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.5)
 * }

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

$$ b''_i=-i^2cos(ix) \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.6)
 * }

Now take the the first and second derivative of Eq. 3.4.2

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

$$ b'_j=-jsin(jx) \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.7)
 * }

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

$$ b''_j=-j^2cos(jx) \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.8)
 * }

Now, substitute values from Eqs. 3.4.1, 3.4.2, and 3.4.5-8
 * {| style="width:100%" border="0"

$$\alpha = cos(ix)[a^'_{2}(-jsin(x)) + a_{2}(-j^2cos(jx))]$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.9)
 * }

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

$$\beta = cos(jx)[a^'_2(-isin(x)) + a_2(-i^2cos(ix))]$$ $$ Now, compare Eq. 3.4.10 to Eq.3.4.10. Since i≠j
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.10)
 * }

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

$$ cos(ix)[a^'_{2}(-jsin(x)) + a_{2}(-j^2cos(jx))] \ne cos(jx)[a^'_2(-isin(x)) + a_2(-i^2cos(ix))]$$ $$ Therefore,
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.4.11)
 * }

α≠β

--Eml5526.s11.team5.smith 18:30, 15 February 2011 (UTC)

= 3.5 Equivalent stiffness of a spring for the bar of thickness t with a centered square hole =

Problem
Show that the equivalent stiffness of a spring aligned in the x direction for the bar of thickness t with a centered square hole shown in the figure is :

where $$ E \quad $$ is the Young's Modulus and t is the width of the bar

Solution
To find the equivalent stiffness of the bar, we have to subdivide the bar into 4 parts and each part behaves like an individual spring with stiffness K1, K2, K3 and K4 respectively. It can be shown from the figure below.

The Equivalent spring - mass system can be shown from the figure below.

From the Hooke's Law,

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

$$\sigma = E \varepsilon$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.2)
 * }

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

$$\frac{F}{A} = E \frac{\delta }{L}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.3 )
 * }

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

$$\frac{F}{\delta} = E \frac{A}{L}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.4 )
 * }

where,

$$\sigma - \quad$$      Stress on the bar (in Pa or psi)

$$E - \quad$$Young's Modulus of the material(in Pa or psi)

$$\varepsilon - \quad$$ Strain experienced by the bar(no units)

$$F - \quad$$Restoring force exerted on the bar (in SI units: N or kgms-2

$$A - \quad$$Area of the cross section of the bar(in sq.m or sq.in)

$$\delta - \quad$$Displacement of the bar(in m or inch)

$$L - \quad$$Length of the bar(in m or inch)

From the definition of Stiffness(K),

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

$$K = \frac{F}{\delta}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

where,

$$K - \quad$$ Stiffness of the bar (N/m of lbf/in)

Therefore $$ \displaystyle (Eq. 3.5.4 ) $$ becomes

$$K = E \frac{A}{L}$$ $$
 *
 $$ \displaystyle (Eq. 3.5.5   )

Hence the stiffness of the springs K1, K2, K3 and K4 can be determined as,

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

$$K_{1} = \frac{EA_{1}}{L_{1}} = \frac{Eat}{\frac{l}{10}}$$                   (From figure 2) $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ K_{1}= \frac{10Eat}{l}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.6  )
 * }

Similarly

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

$$K_{4}= \frac{10Eat}{l}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq.3.5.7   )
 * }

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

$$K_{2} = \frac{EA_{2}}{L_{2}} = \frac{E. 2b. t}{\frac{4l}{5}}$$                  (From figure 2) $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$K_{2}= \frac{5Ebt}{2l}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.8  )
 * }

Similarly

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

$$K_{3}= \frac{5Ebt}{2l}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.9  )
 * }

From Figure 2 and 3, The total stiffness of the bar (K) can be derived as follows,

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

$$\frac{1}{K} = \frac{1}{K_{1}} + \frac{1}{(K_{2} + K_{3})} + \frac{1}{K_{4}}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.10)
 * }

Due to symmetry,$$K_{1} = K_{4}\quad$$ and $$K_{2} = K_{3},\quad$$ therefore,

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

$$\frac{1}{K} = \frac{2}{K_{1}} + \frac{1}{2K_{2}}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.5.11)
 * }

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

$$K = \frac{2K_{1}K_{2}}{K_{1}+4K_{2}} = \frac{2(\frac{10Eat}{l})(\frac{5Ebt}{2l})}{(\frac{10Eat}{l} + \frac{20Ebt}{2l})}\quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$K = \frac{(\frac{10Eat}{l})(\frac{5Ebt}{l})}{\frac{10Eat}{l}(a+b)}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

$$K = \frac{5Etab}{(a+b)l}$$

--Eml5526.s11.team5.srv 02:15, 13 February 2011 (UTC)

= 3.6 : 2D Truss Problem=

Given
a 3-Bar structure is subjected to a load P(=1000N) at C as shown in Figure 3.3.1.

The cross sectional areas are given as:

BC: 0.02

BF,BD: 0.01

All the bars are made of the same material whose Youngs modulus $$ E = 10^{11} $$ Pa.

Find
a) Construct the global stiffness matrix and Load matrix

b) Partition the matrices and solve for the x,y displacements at B, x displacement at D.

c) Find the stress in the 3 bars BC, BD, BF.

d) Find the Reactions at the nodes C, D, F.

Solution
$$ \mathbf {a)} $$

The structure can be divided into 3 elements:

Element 1 = Bar BF

Element 2 = Bar BC

Element 3 = Bar BD (Element numbers are shown in the above figure)

Before we construct the Stiffness matrix for each bar element, let us review the element stiffness matrix for a bar element having 2 degrees of freedom (X,Y displacements) and inclined at an angle $$ \phi $$ with the X-axis.

$$ K^{e}=k^{e} \begin{bmatrix} \cos^2 \phi^e & \cos \phi^e \sin \phi^e & -\cos^2 \phi^e  & -\cos \phi^e \sin \phi^e \\ \cos \phi^e \sin \phi^e & \sin^2 \phi^e & -\cos^2 \phi^e &-\sin^2 \phi^e \\ -\cos^2 \phi^e & -\cos \phi^e \sin \phi^e & \cos^2 \phi^e  & \cos \phi^e \sin \phi^e \\ -\cos \phi^e \sin \phi^e & -\sin^2 \phi^e & \cos^2 \phi^e &\sin^2 \phi^e \\ \end{bmatrix} $$ (ref: Fish and Belytschko P30, Eq (2.47))

where

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

$$ k^e= \frac{A^eE^e}{l^e} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

$$A^e=$$ Cross sectional area of the bar element

$$E^e=$$ Young's Modulus of the bar element

$$l^e=$$ Length of the bar element

$$ \mathbf {Element 1:}$$

Element 1 connects global nodes 1,4.

$$\phi^{(1)}= \frac{\pi}{4}$$

$$\quad A^{(1)}= 0.01$$

$$\quad E^{(1)}=10^{11}$$

$$l^{(1)}=\sqrt{2}$$

$$ \implies K^{(1)}=\frac{10^{9}}{\sqrt{2}} \begin{matrix}

\begin{bmatrix} \frac{1}{2} & \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2}\\ \frac{1}{2} & \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2}\\ -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2}\\ -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \end{bmatrix}

\begin{matrix} \\ (1)\\ \\ (2)\\ \\ \end{matrix}

\\

\begin{matrix} &(1) & & (2)  & \end{matrix}

\end{matrix} $$

$$ \mathbf {Element 2:} $$

Element 2 connects global nodes 1,3.

$$\phi^{(2)}= \frac{\pi}{2}$$

$$\quad A^{(2)}= 0.02$$

$$\quad E^{(2)}=10^{11}$$

$$\quad l^{(2)}=1$$

$$ \implies K^{(2)}=2*10^{9} \begin{matrix}

\begin{bmatrix} 0 & 0 & 0 & 0\\ 0 & 1 & 0 & -1\\ 0 & 0 & 0 & 0\\ 0 & -1 & 0 & 1\\ \end{bmatrix}

\begin{matrix} \\ (1)\\ \\ (3)\\ \\ \end{matrix}

\\

\begin{matrix} &(1) & & (3)  & \end{matrix}

\end{matrix} $$

$$ \mathbf {Element 3:} $$

Element 3 connects global nodes 1,4.

$$\phi^{(3)}= \frac{3\pi}{4}$$

$$\quad A^{(3)}= 0.01$$

$$\quad E^{(3)}=10^{11}$$

$$l^{(3)}=\sqrt{2}$$

$$ \implies K^{(3)}=\frac{10^{9}}{\sqrt{2}} \begin{matrix}

\begin{bmatrix} \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2}\\ -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} & -\frac{1}{2}\\ -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} & -\frac{1}{2}\\ \frac{1}{2} & -\frac{1}{2} & -\frac{1}{2} & \frac{1}{2} \end{bmatrix}

\begin{matrix} \\ (1)\\ \\ (4)\\ \\ \end{matrix}

\\

\begin{matrix} &(1) & & (4)  & \end{matrix}

\end{matrix} $$

Direct Assembly Gives the global matrices:

$$ \implies K={10^{9}} \begin{matrix}

\begin{bmatrix} \frac{1}{\sqrt2}& 0 & -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & 0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2}\\ 0 & 2+\frac{1}{\sqrt2} & -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & -2 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} \\ -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 \\ -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & -2 & 0 & 0 & 0 & 2 & 0 & 0\\ -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} \\ \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \\ \end{bmatrix}

\begin{matrix} \\ (1)\\ \\ (2)\\ \\ (3)\\ \\ (4)\\ \\ \end{matrix}

\\

\begin{matrix} (1)& & &(2)& & &(3)& & &(4) \end{matrix}

\end{matrix} $$

$$ d=\begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{2x}\\ 0\\ 0\\ 0\\ 0\\ 0 \end{bmatrix},

f=1000\begin{bmatrix} 1\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0 \end{bmatrix},

r=\begin{bmatrix} 0\\ 0\\ 0\\ r_{2y}\\ r_{3x}\\ r_{3y}\\ r_{4x}\\ r_{4y} \end{bmatrix} $$

$$ \mathbf {b)} $$

$$ \implies

10^9\begin{bmatrix} \frac{1}{\sqrt2}& 0 & -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & 0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2}\\ 0 & 2+\frac{1}{\sqrt2} & -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & -2 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} \\ -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 \\ -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & -2 & 0 & 0 & 0 & 2 & 0 & 0\\ -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} \\ \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \\ \end{bmatrix}

\begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{2x}\\ 0\\ 0\\ 0\\ 0\\ 0 \end{bmatrix}

=

\begin{bmatrix} 1000\\ 0\\ 0\\ r_{2y}\\ r_{3x}\\ r_{3y}\\ r_{4x}\\ r_{4y} \end{bmatrix} $$

It can be observed that partitioning needs to be done after 3 rows and 3 columns.

The reduced system of equations are:

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

 * style="width:95%" |
 * style="width:95%" |

$$ K_{F}= 10^9\begin{bmatrix} \frac{1}{\sqrt2}& 0 & -\frac{1}{2\sqrt2}\\ 0 & 2+\frac{1}{\sqrt2} & -\frac{1}{2\sqrt2}\\ -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \end{bmatrix},

K_{EF}= 10^9\begin{bmatrix} -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \\ 0 & 0 & 0\\ 0 & -2 & 0\\ -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 \\ \frac{1}{2\sqrt2}& -\frac{1}{2\sqrt2} & 0 \end{bmatrix},

K_{E}= 10^9\begin{bmatrix} \frac{1}{2\sqrt2} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 2 & 0 & 0\\ 0 & 0 & 0 & \frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} \\ 0 & 0 & 0 & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \end{bmatrix}, $$

$$
 *
 $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ d_{F}= \begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{2x} \end{bmatrix},

\bar d_{E}= \begin{bmatrix} 0\\ 0\\ 0\\ 0\\ 0 \end{bmatrix} = \begin{bmatrix} \bar u_{2y}\\ \bar u_{3x}\\ \bar u_{3y}\\ \bar u_{4x}\\ \bar u_{4y} \end{bmatrix},

f_{F}= \begin{bmatrix} 1000\\ 0\\ 0 \end{bmatrix},

r_{E}= \begin{bmatrix} r_{2y}\\ r_{3x}\\ r_{3y}\\ r_{4x}\\ r_{4y} \end{bmatrix}, $$

$$
 *
 $$ \displaystyle\!
 * }

$$ \implies \begin{bmatrix} K_{F} & K_{EF}^{T}\\ K_{EF} & K_{E} \end{bmatrix}

\begin{bmatrix} d_{F}\\ \bar d_{E} \end{bmatrix} = \begin{bmatrix} f_{F}\\ r_{E} \end{bmatrix} $$

In other words:

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

$$K_{F} d_{F} + K_{EF}^T \bar d_{E} = f_{F}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq 3.6.1)
 * }

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

$$K_{EF} d_{F} + K_{E} \bar d_{E} = r_{E}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq 3.6.2)
 * }

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

$$d_{E}=\begin{bmatrix} 0&0&0&0&0\end{bmatrix}^{T}$$, $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq 3.6.3)
 * }

Combining Eq 3.6.1, Eq 3.6.3:
 * {| style="width:100%" border="0"

$$\quad K_{F} d_{F} = f_{F}$$ $$
 * style="width:100%" |
 * style="width:100%" |
 *
 $$ \displaystyle (Eq 3.6.4)
 * }

Combining Eq 3.6.2, Eq 3.6.3:
 * {| style="width:100%" border="0"

$$\quad K_{EF} d_{F} = r_{E}$$ $$
 * style="width:100%" |
 * style="width:100%" |
 *
 $$ \displaystyle (Eq 3.6.5)
 * }

from Eqn 3.6.4:

$$d_{F} = K_{F}^{-1} f_{F}$$

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

$$ \therefore d_{F}= \begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{2x} \end{bmatrix} = \frac{1}{10^6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 1+2\sqrt{2} \end{bmatrix} $$

X,Y-displacement at B(in mm) $$\quad = \left( u_{1x}, u_{1y} \right)=\left( 3.3284*10^{-3}, 5*10^{-4} \right) $$

X-displacement at D(in mm) $$\quad = \left( u_{2x} \right)=\left( 3.8284*10^{-3} \right) $$

$$ \mathbf {c)} $$

The stress in each element can be computed as follows:

$$ \sigma^e = E^e \frac{u_{2x}^{'e}-u_{1x}^{'e}}{l^{e}} = \frac{E^e}{l^e} \begin{bmatrix} -1 & 0 & 1 & 0 \end{bmatrix} \begin{bmatrix} u_{1x}^{'e}\\ u_{1y}^{'e}\\ u_{2x}^{'e}\\ u_{2y}^{'e} \end{bmatrix} =\frac{E^e}{l^e} \begin{bmatrix} -\cos \phi^e & -\sin \phi^e & \cos \phi^e & \sin \phi^e \end{bmatrix} \mathbf d^e $$ (ref: Fish and Belytschko P30, Eq (2.47))

Computing the Stresses in each bar (element) using the above formula:

$$ \mathbf {Element 1:}$$

$$ \mathbf d^{(1)}= \begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{2x}\\ u_{2y} \end{bmatrix} = \frac{1}{10^6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 1+2\sqrt{2}\\ 0 \end{bmatrix} $$

$$ \implies \sigma^{(1)}=\frac{10^{11}}{\sqrt{2}} \begin{bmatrix} -\frac{1}{\sqrt 2} & -\frac{1}{\sqrt 2} & \frac{1}{\sqrt 2} & \frac{1}{\sqrt 2} \end{bmatrix}* \frac{1}{10^6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 1+2\sqrt{2}\\ 0 \end{bmatrix} $$

$$ \implies \sigma^{(1)}=0$$

$$ \mathbf {Element 2:}$$

$$ \mathbf d^{(2)}= \begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{3x}\\ u_{3y} \end{bmatrix} = \frac{1}{10^6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 0\\ 0 \end{bmatrix} $$

$$ \implies \sigma^{(2)}=10^{11} \begin{bmatrix} 0 & -1 & 0 & 1 \end{bmatrix}* \frac{1}{10^6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 0\\ 0 \end{bmatrix}$$

$$ \implies \sigma^{(2)}=-50*10^3 = -50KPa$$

$$ \mathbf {Element 3:}$$

$$ \mathbf d^{(3)}= \begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{4x}\\ u_{4y} \end{bmatrix} = \frac{1}{10^6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 0\\ 0 \end{bmatrix} $$

$$ \implies \sigma^{(3)}=\frac{10^{11}}{\sqrt{2}} \begin{bmatrix} \frac{1}{\sqrt 2} & -\frac{1}{\sqrt 2} & -\frac{1}{\sqrt 2} & \frac{1}{\sqrt 2} \end{bmatrix}* \frac{1}{10^6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 0\\ 0 \end{bmatrix} $$

$$ \implies \sigma^{(3)}=\sqrt{2} * 10^5Pa= 141.42KPa$$

$$ \mathbf {d)} $$

Substituting $$ d_{F} $$ in Eqn 3.6.2, we calculate the reactions forces ($$r_{E}$$ matrix)

$$

r_{E}= \begin{bmatrix} r_{2y}\\ r_{3x}\\ r_{3y}\\ r_{4x}\\ r_{4y} \end{bmatrix} = 10^9\begin{bmatrix} -\frac{1}{2\sqrt2} & -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} \\ 0 & 0 & 0\\ 0 & -2 & 0\\ -\frac{1}{2\sqrt2} & \frac{1}{2\sqrt2} & 0 \\ \frac{1}{2\sqrt2}& -\frac{1}{2\sqrt2} & 0 \end{bmatrix} \frac{1}{10^{6}} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 1+2\sqrt{2} \end{bmatrix} $$

$$ r_{E}= \begin{bmatrix} r_{2y}\\ r_{3x}\\ r_{3y}\\ r_{4x}\\ r_{4y} \end{bmatrix} =10^3 \begin{bmatrix} 0\\ 0\\ -1\\ -1\\ 1 \end{bmatrix} $$

The reactions at each Node (in N):

C(X-component,Y-component) = (0,-1000)

D(Y-component) = 0

F(X-component,Y-component) = (-1000, 1000)

Comparing Results with ABAQUS
The Stresses and Displacements calculated have been compared with the results generated by ABAQUS and there is a very good match as seen below.

2D Truss Model in ABAQUS, after the Loads and BCs have been applied.

The above image displays the Stress in each bar.

The above image displays the X-Displacements of the structure.

The above image displays the Y-Displacements of the structure.

--Eml5526.s11.team5.vijay 19:56, 15 February 2011 (UTC)

= 3.7 Weighted Residue Form as BC on Polynomial Base Function (k = 0) =

Given

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

$$ \omega = ]0,1[ \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.1)
 * }

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

$$ a_2 = 2 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.2)
 * }

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

$$ f = 3 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.3)
 * }

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

$$ \frac{\partial{u}}{\partial{t}}(x,t) = 0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.4)
 * }

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

$$ \Gamma_g = 1 $$, $$ g = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ u(x=1) = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.5)
 * }

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

$$ \Gamma_h = 0 $$, $$ h = 4 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ -\frac{du^h}{dx}(0) = 4 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.7.6)
 * }

Find
 Find approximate solution $$ u^n \quad $$ and compare it to exact one

Consider base function $$ b_j(x), j = 0, 1, 2, 3,... n \quad $$

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

$$ b_j(x) = (x+k)^j \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.7)
 * }

,where $$ k \quad $$ is an arbitrary constant, which in this exercise it chosen to be $$ 0 \quad $$

(1) Let $$ n = 2 \Rightarrow ndof = n+1 = 3 \quad $$

Generate variables considered for $$ n = 2 \quad $$ case

(2) Find 2 eqs that enforce b.c.'s for $$ u^h(x) = \sum_{j -0}^{n} d_j b_j(x) $$

(3) Find 1 more eq. to solve for $$ \underline{d} = d_j (j = 0, 1, 2) $$ by projecting the residue $$ \underline{p}(u^h) $$ on a basis function $$ b_k(x) \quad $$ with $$ k = 0, 1, 2 \quad $$. The additional equation must be lineally independent from the above 2 equations in part (2).

(4) Display 3 equations in matrix form

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

$$ \underline{K}\underline{d} = \underline{F} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.8)
 * }

Observe if $$ \underline{K} $$ is symmetric

(5) Solve for $$ \underline{d} \quad $$

(6) Construct $$ u_n^k(x) $$ and plot $$ u_n^k(x) $$ vs. $$ ux) $$

(7) Repeat steps (1) to (6) for

 * (7.1) $$ n = 4 \quad $$
 * (7.2) $$ n = 6 \quad $$

(8) Compute $$ u_n^h (x = 0.5) \quad $$ for $$ n = 2, 4, 6 \quad $$

error

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

$$ e_n(0.5) = u_n(0.5) - u_n^h(0.5) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.9)
 * }

Plot $$ e_n(0.5) \quad $$ vs. $$ n \quad $$

Solution
(1)

For $$ n = 2 \quad $$

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

$$ \underline{b} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} b_0     & b_1 & b_2 \end{bmatrix}

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

From base function expression $$ Eq. 3.7.7 \quad $$

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

$$ \underline{b} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} (x+0)^0     & (x+0)^1 & (x+0)^2 \end{bmatrix}

=

\begin{bmatrix} 1     & x & x^2 \end{bmatrix}

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

Unknown function

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

$$ \underline{d} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} d_0     & d_1 & d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle (Eq. 3.7.12)
 * }

Approximate solution

(2)

at B.C. $$ u(1) = 0 \quad $$

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

$$ u^h = \sum_{i=0}^{n=2} d_ib_i $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.13)
 * }

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

$$ u^h(x) = d_0(1)   + d_1x +d_2x^2   \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.14)
 * }

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

$$ u^h(1) = 0 = d_0(1)   + d_1(1) +d_2(1)^2 =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1     & 1 & 1      \end{bmatrix}

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle\!
 * }

at B.C $$ \frac{du}{dx}(0) = 4 $$

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

$$ (u')^h(x) = \frac{du}{dx} = \sum_{i=0}^{n=2} d_ib_i' $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.15)
 * }

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

$$ \underline{b}' =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0     & 1 & 2x \end{bmatrix}

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

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

$$ (u')^h(0) = - 4 = d_00 + d_11 + d_22(0) =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0     & 1 & 0    \end{bmatrix}

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle\!
 * }

(3)

By projecting the residue part of a function on the bases, 3rd equation can be generated:

From the definition of elastodynamics case of ODE function $$ P_u(x) \quad $$

Rewriting $$ Eq. 3.7.7 \quad $$

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

$$ b_j(x) = (x+k)^j \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.7)
 * }

In this, more specific, "self adjoin" case

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

$$ P(u^h) :=[a_2u']' + a_0u - g = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.17)
 * }

noting that

$$ g = 0 \Rightarrow \quad $$ no forcing function (steady state)

$$ a_0u = f = 3 \quad $$ (given)

$$ a_2 = 2 \quad $$ (given)

$$ P(u^h) :=[2u']' + 3 = 0 \Rightarrow \quad $$ since $$ a_2 \quad $$ is given as constant 2 (not function of "x"), it can be taken out of bracketing differential, resulting in

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

 * style="width:95%" |
 * style="width:95%" |

$$ P(u^h) :=2u'' + 3 = 0= 2\frac{d^2u}{dx} +3 = 0 $$ $$
 *
 $$ \displaystyle (Eq. 3.7.18)
 * }

From $$ Eq. 3.7.13 \quad $$, $$ u^h \quad $$ was defined as

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

$$ u^h = \sum_{i=0}^{n=2} d_ib_i = u^h(x) = d_0(1)   + d_1x +d_2x^2   \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.19)
 * }

Taking $$ 2^{nd} \quad $$ derivative with respect to "x"

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

$$ u^{h''} = d_00 + d_10 + 2d_2 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.20)
 * }

Therefore

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

$$ P(u^h) := 2(d_00 + d_10 + 2d_2)) +3 = 4d_2 + 3\quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.21)
 * }

Afterwords, noting that

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

$$ \int_\omega w P(u^h) dx := 0 \Rightarrow w_i = b_i \Rightarrow \int_\omega b_i P(u^h) dx := 0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.22)
 * }

which is also noted in $$ 10-4, Eq. (2) \quad $$

and substituting for $$ P(u^h) \quad $$ equation above, we get

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

 * style="width:95%" |
 * style="width:95%" |

$$ \int_\omega b_i (4d_2 + 3) dx := 0 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.7.23)
 * }

recalling that $$ \underline{b} =

\begin{bmatrix} 1  \\ x \\ x^2 \end{bmatrix} $$

and taking any of the equations as the last needed equation, for example $$ b_1 \quad $$, we get

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

$$ \int_\omega
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1  \\ x \\ x^2 \end{bmatrix}

[0d_0 + 0d_1 + 4d_2] dx := \int_\omega (-3

\begin{bmatrix} 1  \\ x \\ x^2 \end{bmatrix} )dx \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.7.25)
 * }

Doing symbolic integration and evaluating interval on $$ \omega = ]0,1[ \quad $$ we get

The outputs are as follow:

$$ \int_\omega

\begin{bmatrix} 1  \\ x \\ x^2 \end{bmatrix}

[0d_0 + 0d_1 + 4d_2] dx = $$ 0           0  4.0000e+000 0           0  2.0000e+000 0           0  1.3333e+000

$$ = \int_\omega (-3

\begin{bmatrix} 1  \\ x \\ x^2 \end{bmatrix} )dx \quad $$

-3.0000e+000 -1.5000e+000 -1.0000e+000

(4)

Assembling the matrix $$ \underline{K} \underline{d} = \underline{F} \quad $$ from 2 B.C. constrains and 1 more equation generated by weighting factor method:

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

$$ \begin{bmatrix} 1     & x & x^2 \\ 0 &   1  & 2x \\ 0 &   0  & 4     \end{bmatrix}
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix}

=

\begin{bmatrix} 0\\ 4 \\ -3   \end{bmatrix} $$ $$
 *
 $$ \displaystyle (Eq. 3.7.26)
 * }

For "K" matrix: 1    1     1     0     1     0     0     0     4

For "F" matrix:

0   -4    -3

(5)

Solving for $$ \underline{d} = \underline{K}^{-1} F\quad $$ we get

d = 4.7500e+000 -4.0000e+000 -7.5000e-001

(6)

Constructed $$ u^h \quad $$ matrix on matlab by using formula $$ u^h = \sum_{i=0}^{n=2} d_ib_i $$

4.75 - 4.0*x - 0.75*(k + x)^2 - 4.0*k

Exact solution:

Given a general expression for 2nd Order ODE

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

$$ a_2(x)y'' + a_1(x)y + a_0(x)y = g(x) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.27)
 * }

and taking a more specific case ("Self adjoin") for this specific problem, we have

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

$$ [a_2(x)y']' + a_0(x)y = g(x) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.17)
 * }

$$ g(x) = 0 \Rightarrow \quad $$ no forcing function (steady state)

$$ a_0(x)y = f(x) = 3 \quad $$ (given)

$$ a_2 = 2 \quad $$ (given)

Therefore $$ Eq. 3.7.17 \quad $$ becomes

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

 * style="width:95%" |
 * style="width:95%" |

$$ [2y']' + 3 = 0 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.7.28)
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ \frac{du}{dx}(\frac{du}{dx}2) + 3 = 0 $$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ \int \frac{du}{dx}2 = \int -3 dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \frac{du}{dx}2 = -3x + c_1' \Rightarrow \frac{du}{dx} = -1.5x + c_1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.29)
 * }

at B.C of $$ -\frac{du}{dx}(0) = 4 \Rightarrow $$

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

 * style="width:95%" |
 * style="width:95%" |

$$ -4 = -1.5(0) + c_1 \Rightarrow c_1 = -4$$ $$
 *
 $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ \int {du} = \int (-1.5x - 4) dx $$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ u = \frac{-1.5}{2}x^2 - 4x + c_2 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.30)
 * }

at B.C of $$ u(1) = 0 \Rightarrow $$

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

 * style="width:95%" |
 * style="width:95%" |

$$ 0 = \frac{-1.5}{2}(1)^2 - 4(1) + c_2 \Rightarrow c_2 = -4.75$$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ u(x) = -0.75x^2 - 4x - 4.75 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.31)
 * }

Circle view

So taken in account the exact solution formula, different plot were generated for different n values (n = 2, 4, 6)

Error came out to be:

x = 0.5

ex = n=2       n=4          n=6 0           0      4.3521e-014

Note, that the 2 member bases approximated the given problem more exact then n = 4 or n = 6. To be exact, it approximated exactly. This was unusual, it might have had been a case that the point around which the error was calculated just happen to lay at that point. So we checked for x = 1, x = 0.25 and x = 0.75. Here are the results:

x = 1

ex = n=2       n=4          n=6 0          0       1.0800e-012

x = 0.25

ex = n=2       n=4          n=6 0        0       1.3323e-015

x = 0.75

ex = n=2       n=4          n=6 0         0       2.8599e-013

What results basically advocated was that n =2,4 exactly approximated 2nd order ODE, while for this base, increasing number of members 'n' produce relatively greater error (see graphs below). Also interesting observation could be made that with k = 0 value, n = 2 and n = 4 exactly approximates the polynomial while with k $$ \ne $$ to 0, only n=2 can exactly approximate the ODE (see exercise 1 for k $$ \ne $$ to 0 case).

Netherless, from previous experiences such as convergence of solutions that are seen below and re-use of the earlier proven to work code, it is clear that for these bases $$ b_j = (x+k)^j \quad $$ higher order terms induce relatively higher error then lower order terms.

Approximation and error outputs for bases $$ b_j(x) = (x+k)^j \quad $$

Approximation zoomed in view 251 equidistant point were taken for non-symbolic plot of n=2,4,6, where black single dot represent numerical n=2 value at x=0.5 solution and red circle represents symbolic value (blue line) of exact solution at x = 0.5

--Egm6341.s11.team3.JA 02:32, 10 February 2011 (UTC)

=3.8 Derivation of Weak Form from Strong Form =

Given

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

The strong form is given by: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \frac{d}{dx} \left (AE \frac{du}{dx} \right ) +2x=0 on 1 < x < 3 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.1)
 * }

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

$$ \sigma (1)= \left (E \frac{du}{dx} \right )_{x=1}=0.1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.2)
 * }

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

$$ u(3)=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.3)
 * }

Show

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

The weak form is given by: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \int\limits_{1}^{3}\frac{dw}{dx} AE \frac{du}{dx} dx= -0.1(wA)_{x=1}+ \int\limits_{1}^{3} 2xw dx ~\forall w ~  with ~ w(3)=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.4)
 * }

Solution
Multiply both Eq. 3.8.1 and 3.8.2 by arbitrary function w(x) and integrate over their domains. (Note: Eq. 3.8.2 does not need an integral since the condition only exists at one point; instead, multiply by area 'A')

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

$$\int\limits_{1}^{3} w [\frac{d}{dx} \left (AE \frac{du}{dx} \right ) +2x]dx=0 ~\forall w$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.5)
 * }

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

$$(wA(E \frac{du}{dx} -.1))_{x=1}=0~\forall w$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.6)
 * }

The w(x) is the weight function, and it is convenient to have it satisfy

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

$$ w(3)=0 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.6)
 * }

Eq. 3.8.5 then becomes

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

$$\int\limits_{1}^{3} w \frac{d}{dx} \left (AE \frac{du}{dx} \right )dx +\int\limits_{1}^{3} w2xdx=0 ~\forall w$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.7)
 * }

Integrate the left side of Eq. 3.8.7, by parts

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

$$\int\limits_{1}^{3} w \frac{d}{dx} \left (AE \frac{du}{dx} \right )dx = (wAE \frac{du}{dx})_{x=1,3}- \int\limits_{1}^{3}  \frac{dw}{dx}AE \frac{du}{dx} dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

Now write Eq. 3.8.7

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

$$ (wA\underbrace{E \frac{du}{dx}}_\sigma)_{x=1,3} - \int\limits_{1}^{3} \frac{dw}{dx}AE \frac{du}{dx} dx - \int\limits_{1}^{3} w2xdx=0 ~\forall w, w(3)=0$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.8)
 * }

By the stress-strain law (See Hooke's law)and strain displacement (See Strain)equations, the above eqn can be simplified
 * {| style="width:100%" border="0"

$$(wA\sigma)_{x=3} = 0 ; (wA\sigma)_{x=1} = 0.1$$ $$ Eq. 3.8.8 then becomes
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }
 * {| style="width:100%" border="0"

$$0-.1- \int\limits_{1}^{3} \frac{dw}{dx}AE \frac{du}{dx} dx - \int\limits_{1}^{3} w2xdx=0 ~\forall w, w(3)=0$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.9)
 * }

$$ \int\limits_{1}^{3}\frac{dw}{dx} AE \frac{du}{dx} dx= -0.1(wA)_{x=1}+ \int\limits_{1}^{3} 2xw dx ~\forall w ~  with ~ w(3)=0 $$

--Eml5526.s11.team5.smith 18:30, 15 February 2011 (UTC)

= 3.9 Linear trial solution of the weak form =

Given

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

The strong form is given by: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \frac{d}{dx} \left (AE \frac{du}{dx} \right ) +2x=0 on 1 < x < 3 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.9.1)
 * }

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

$$ \sigma (1)= \left (E \frac{du}{dx} \right )_{x=1}=0.1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.9.2)
 * }

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

$$ u(3)=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.9.3)
 * }

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

The weak form of this strong form is given by: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \int\limits_{1}^{3}\frac{dw}{dx} AE \frac{du}{dx} dx= -0.1(wA)_{x=1}+ \int\limits_{1}^{3} 2xw dx ~\forall w ~  with ~ w(3)=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.9.4)
 * }

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

Consider a trial solution of the form $$u(x)= \alpha_0+\alpha_1(x-3) \quad $$ and a weight function of the same form. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

Find

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

1) Obtain a solution to the weak form given. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

2) Check that the equilibrium equation in the strong form is satisfied. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

3) Check that the natural boundary condition is satisfied. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

Solution

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

1) Obtain a solution to the weak form given. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

To be admissible the weight function must vanish at x=3: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$w(x=3)=0 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$w(x=3)=\Beta_0+\Beta_1(3-3)=0 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

This results in $$\Beta_0=0\quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

For the trial solution to be admissible, it must satisfy the essential boundary condition $$u(3)=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$u(x=3)=\alpha_0+\alpha_1(3-3)=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

This gives the result $$\alpha_0=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Substituting these results back into the weight function and the trial solution yields the following equations and their derivatives: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ u(x)=0.001+\alpha_1(x-3),~\frac{du(x)}{dx}=\alpha_1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.9.5)
 * }

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

$$ w(x)=\Beta_1(x-3),~\frac{dw(x)}{dx}=\Beta_1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.9.6)
 * }

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

Substituting Eq. 3.9.5 and Eq. 3.9.6 into the weak form (Eq. 3.9.4) gives: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \int\limits_{1}^{3}\Beta_1 \alpha_1 AE dx =-0.1(\Beta_1(x-3)A)_{x=1} + \int\limits_{1}^{3}2x \Beta_1 (x-3)dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Assuming the area, A and young's modulus, E are constants yields: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \Beta_1 AE \alpha_1 \int\limits_{1}^{3}dx = -0.1(\Beta_1 (x-3)A)_{x=1}+2 \Beta_1 \int\limits_{1}^{3} x(x-3) dx$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Evaluating the integrals yields: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \Beta_1 \alpha_1 AE (3-1)=\frac{1}{5} \Beta_1 A +2 \Beta_1 \left ( -\frac{9}{2}+\frac{7}{6}\right )$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Simplifying and moving everything to the left gives: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ 2\Beta_1\alpha_1AE= \frac{1}{5} \Beta_1 A -\frac{20}{3} \Beta_1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ 2\Beta_1\alpha_1 AE-\frac{1}{5} \Beta_1 A +\frac{20}{3} \Beta_1 =0$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \Beta_1 \left (2\alpha_1 AE-\frac{1}{5}A +\frac{20}{3} \right )=0$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.9.7)
 * }

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

Since Eq. 3.9.7 most hold for all $$ \Beta_1$$ the term in the parentheses must be equal to 0: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ 2 \alpha_1 AE-\frac{1}{5}A +\frac{20}{3} =0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Solving for $$ \alpha_1$$ gives: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\alpha_1=\frac{3A-100}{30AE}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Substituting this result into Eq. 3.9.5 gives the weak form solution: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$u(x)=0.001+\frac{3A-100}{30AE}(x-3)$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

2) Check the equilibrium equation in the strong form $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

In order to check the equilibrium equation in the strong form Eq. 3.9.1 is used with the values we just found: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\frac{d}{dx} \left (AE \frac{3A-100}{30AE} \right )+2x=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\frac{d}{dx} \left (\frac{3A-100}{30}\right )+2x=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$0+2x=0 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

 :::{| style="width:100%" border="0" As we can see the equilibrium in the strong form is not satisfied. Thus this trial solution will not work and another form will need to be used. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

3) Check that the natural boundary condition is satisfied. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

In order to check the natural boundary condition Eq. 3.9.2 is used: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \left (E \frac{du}{dx} \right )_{x=1}=0.1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \left (E \frac{3A-100}{30AE} \right )_{x=1}=0.1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \left (\frac{3A-100}{30A} \right )=0.1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ 3A-100 =3A \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

We can clearly see that this is not possible and that the natural boundary condition is not met with this type of trial solution. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

= 3.10 Quadratic trial solution of the weak form =

Given

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

The strong form is given by:
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!$$
 * }
 * {| style="width:100%" border="0"

$$ \frac{d}{dx} \left (AE \frac{du}{dx} \right ) +2x=0 on 1 < x < 3 $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.10.1)$$
 * }

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

$$ \sigma (1)= \left (E \frac{du}{dx} \right )_{x=1}=0.1 $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.10.2)$$
 * }

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

$$ u(3)=0.001 \quad$$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.10.3)$$
 * }

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

The weak form of this strong form is given by:
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!$$
 * }

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

$$ \int\limits_{1}^{3}\frac{dw}{dx} AE \frac{du}{dx} dx= -0.1(wA)_{x=1}+ \int\limits_{1}^{3} 2xw dx ~\forall w ~  with ~ w(3)=0 $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.10.4)$$
 * }

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

Consider a trial solution of the form $$u(x)= \alpha_0+\alpha_1(x-3)+\alpha_2(x-3)^{2} \quad $$ and a weight function of the same form.
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!$$
 * }

Find

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

1) Obtain a solution to the weak form given. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

2) Check that the equilibrium equation in the strong form is satisfied. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

3) Check that the natural boundary condition is satisfied. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

Solution

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

1) Obtain a solution to the weak form given. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

To be admissible the weight function must vanish at x=3: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$w(x=3)=0 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$w(x=3)=\Beta_0+\Beta_1(3-3)+\Beta_2(3-3)^{2}=0 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

This results in $$\Beta_0=0\quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

For the trial solution to be admissible, it must satisfy the essential boundary condition $$u(3)=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$u(x=3)=\alpha_0+\alpha_1(3-3)+\alpha_2(3-3)^{2}=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

This gives the result $$\alpha_0=0.001 \quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Substituting these results back into the weight function and the trial solution yields the following equations and their derivatives: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ u(x)=0.001+\alpha_1(x-3)+\alpha_2(x-3)^{2},~\frac{du(x)}{dx}=\alpha_1+\alpha_2(2x-6) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.10.5)
 * }

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

$$ w(x)=\Beta_1(x-3)+\Beta_2(x-3)^{2},~\frac{dw(x)}{dx}=\Beta_1+\Beta_2(2x-6) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.10.6)
 * }

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

Substituting Eq. 3.10.5 and Eq. 3.10.6 into the weak form (Eq. 3.10.4) gives: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \int_{1}^{3}(B_1+B_2(2x-6))AE(\alpha_1+\alpha_2(2x-6))dx=-0.1(A)(B_1(x-3)+B_2(x-3)^{2})_{x=1}+ \int_{1}^{3}2x(B_1(x-3)+B_2(x-3)^{2})dx ~\forall w ~  with ~ w(3)=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Now assuming that the area, A and young's modulus, E are constants, the integrals on the LHS and RHS are both expanded, integrated and evaluated within the respective limits of integration. The resulting equation from above is simply reduced to the following form: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ AE[2B_1\alpha_1-B_1\alpha_2(8)-B_2\alpha_1(8)+B_2\alpha_2(\frac{176}{3})]=-0.1(A)(B_1(-2)+B_2(4))+(\frac{-20}{3})B_1+8B_2 ~\forall w ~  with ~ w(3)=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

In factoring out the B terms first on the LHS and then on the RHS yields: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ 2(AE)B_1(\alpha_1-4\alpha_2)+8(AE)B_2(\frac{22}{3}\alpha_2-\alpha_1)=B_1(\frac{A}{5}-\frac{20}{3})+B_2(8-\frac{2A}{5}) ~\forall w ~  with ~ w(3)=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

In moving all terms now to the LHS and setting the RHS to zero, the B terms may then be grouped together again, giving the equation: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$B_1[E(A\alpha_1-4A\alpha_2)-(\frac{3A-100}{30})]+B_2[(-4\alpha_1A+\frac{88}{3}A\alpha_2)E-(4-\frac{A}{5})]=0 ~\forall w ~  with ~ w(3)=0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Since the equation above must hold for arbitrary weight functions, it must hold for arbitrary $$ B_1\quad$$ and $$ B_2\quad$$. Therefore according to the scalar product theorem, $$ B_1\quad$$ and $$ B_2\quad$$ must vanish, yielding the following linear algebraic equation in terms of $$ \alpha_1\quad$$ and $$ \alpha_2\quad$$. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ E\begin{bmatrix} A&-4A \\ -4A& \frac{88}{3}A \end{bmatrix}\begin{bmatrix} \alpha_1\\ \alpha_2 \end{bmatrix}=\begin{bmatrix} \frac{3A-100}{30}\\ \frac{20-A}{5} \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

If we then factor out the $$ A\quad$$ from the first matrix on the LHS and then take the inverse of that matrix, the following linear algebraic matrix equation results to solve for $$ \alpha_1\quad$$ and $$ \alpha_2\quad$$: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \begin{bmatrix} \alpha_1\\ \alpha_2 \end{bmatrix}=\frac{1}{EA(40)} \begin{bmatrix} 88 &12 \\ 12 &3 \end{bmatrix} \begin{bmatrix} \frac{3A-100}{30}\\ \frac{20-A}{5} \end{bmatrix}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Solving now for $$ \alpha_1\quad$$ and $$ \alpha_2\quad$$ and simplifying gives: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\alpha_1=\frac{4}{EA(75)}(3x-115)$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\alpha_2=\frac{1}{EA(200)}(3x-140)$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Substituting these results for $$\alpha\quad$$ into Eq. 3.10.5 gives the weak form solution: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$u(x)=0.001+\frac{4}{EA(75)}(3x-115)(x-3)+\frac{1}{EA(200)}(3x-140)(x-3)^{2}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.10.7)
 * }

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

2) Check the equilibrium equation in the strong form $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

In order to check the equilibrium equation in the strong form Eq. 3.10.1 is used with the derivative of Equation 3.10.7--which is seen below: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\frac{du(x)}{dx}=\frac{4}{EA(75)}[6x-124]+\frac{1}{EA(200)}[9x^{2}-316x+867]$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Solving now Equation 3.10.1, with $$\frac{du(x)}{dx}\quad$$ yields: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\frac{4}{75}[6]+\frac{1}{200}[18x-316]+2x\overset{?}{=}0 on 1 < x < 3 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

Solving now for $$x\quad$$ it can be seen that $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$2.09x=1.26\quad$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\therefore x=\frac{126}{209}\approx 0.60287 which\ is\ not\ on 1 < x < 3 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

 :::{| style="width:100%" border="0" As we can see the equilibrium in the strong form is not satisfied. Thus this trial solution will not work and another form will need to be used. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

3) Check that the natural boundary condition is satisfied. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

In order to check the natural boundary condition Eq. 3.10.2 is used with $$\frac{du(x)}{dx}\quad$$ yielding: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\sigma (1)=\left [\frac{4}{A(75)}(6x-124)+\frac{1}{A(200)}(9x^{2}-316x+867)\right ]_{x=1}\overset{?}{=}0.1$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

If A=1 in this case, $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$\sigma (1)=\frac{-3649}{300}\approx -12.163\neq 0.1$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

On the other hand if we were to solve for $$A\quad$$, it would have to be equal to -121.6, which for an area is not physically possible. Thus, the natural boundary condition is also not met with this type of assumed solution. $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

--Eml5526.s11.team5.srv 02:15, 13 February 2011 (UTC)

= 3.11 Weighted Residue Form as BC on Fourier Base Function =

Given

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

$$ \omega = ]0,1[ \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.1)
 * }

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

$$ a_2 = 2 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.2)
 * }

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

$$ f = 3 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.3)
 * }

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

$$ \frac{\partial{u}}{\partial{t}}(x,t) = 0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.4)
 * }

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

$$ \Gamma_g = 1 $$, $$ g = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle \!
 * }

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

$$ u(x=1) = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.5)
 * }

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

$$ \Gamma_h = 0 $$, $$ h = 4 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ -\frac{du^h}{dx}(0) = 4 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.11.6)
 * }

Find
 Find approximate solution $$ u^n \quad $$ and compare it to exact one

Consider base function $$ b_j(x), j = 0, 1, 2, 3,... n \quad $$

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

$$ b_j(x) = 1 + Fourier \ series = (1, cos(x), sin(x), cos(2x), sin(2x),.....) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.7)
 * }

(1) Let $$ n = 2 \Rightarrow ndof = n+1 = 3 \quad $$

Generate variables considered for $$ n = 2 \quad $$ case

(2) Find 2 eqs that enforce b.c.'s for $$ u^h(x) = \sum_{j -0}^{n} d_j b_j(x) $$

(3) Find 1 more eq. to solve for $$ \underline{d} = d_j (j = 0, 1, 2) $$ by projecting the residue $$ \underline{p}(u^h) $$ on a basis function $$ b_k(x) \quad $$ with $$ k = 0, 1, 2 \quad $$. The additional equation must be lineally independent from the above 2 equations in part (2).

(4) Display 3 equations in matrix form

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

$$ \underline{K}\underline{d} = \underline{F} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.8)
 * }

Observe if $$ \underline{K} $$ is symmetric

(5) Solve for $$ \underline{d} \quad $$

(6) Construct $$ u_n^k(x) $$ and plot $$ u_n^k(x) $$ vs. $$ ux) $$

(7) Repeat steps (1) to (6) for

 * (7.1) $$ n = 4 \quad $$
 * (7.2) $$ n = 6 \quad $$

(8) Compute $$ u_n^h (x = 0.5) \quad $$ for $$ n = 2, 4, 6 \quad $$

error

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

$$ e_n(0.5) = u_n(0.5) - u_n^h(0.5) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.9)
 * }

Plot $$ e_n(0.5) \quad $$ vs. $$ n \quad $$

Solution
(1)

For $$ n = 2 \quad $$

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

$$ \underline{b} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} b_0     & b_1 & b_2 \end{bmatrix}

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

From base function expression $$ Eq. 3.11.7 \quad $$

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

$$ \underline{b} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1     & cos(x) & sin(x) \end{bmatrix}

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

Unknown function

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

$$ \underline{d} =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} d_0     & d_1 & d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle (Eq. 3.11.12)
 * }

Approximate solution

(2)

at B.C. $$ u(1) = 0 \quad $$

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

$$ u^h = \sum_{i=0}^{n=2} d_ib_i $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.13)
 * }

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

$$ u^h(x) = d_0(1)   + d_1cos(x) +d_2sin(x)  \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.14)
 * }

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

$$ u^h(1) = 0 = d_0(1)   + d_1cos(1) +d_2sin(1) =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1   + cos(1) +sin(1) \end{bmatrix}

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle\!
 * }

at B.C $$ \frac{du}{dx}(0) = 4 $$

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

$$ (u')^h(x) = \frac{du}{dx} = \sum_{i=0}^{n=2} d_ib_i' $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.15)
 * }

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

$$ \underline{b}' =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0     & -sin(x) & cos(x) \end{bmatrix}

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

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

$$ (u')^h(0) = - 4 = d_00 + d_1(-sin(x)) + d_2cos(x) =
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 0     & 0 & 1     \end{bmatrix}

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix} $$ $$
 *
 $$ \displaystyle\!
 * }

(3)

By projecting the residue part of a function on the bases, 3rd equation can be generated:

From the definition of elastodynamics case of ODE function $$ P_u(x) \quad $$

Rewriting $$ Eq. 3.11.7 \quad $$

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

$$ b_j(x) = [1 \ cos(x) \ sin(x)] \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.7)
 * }

for j = 0,1,2

In this, more specific, "self adjoin" case

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

$$ P(u^h) :=[a_2u']' + a_0u - g = 0 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.17)
 * }

noting that

$$ g = 0 \Rightarrow \quad $$ no forcing function (steady state)

$$ a_0u = f = 3 \quad $$ (given)

$$ a_2 = 2 \quad $$ (given)

$$ P(u^h) :=[2u']' + 3 = 0 \Rightarrow \quad $$ since $$ a_2 \quad $$ is given as constant 2 (not function of "x"), it can be taken out of bracketing differential, resulting in

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

 * style="width:95%" |
 * style="width:95%" |

$$ P(u^h) :=2u'' + 3 = 0= 2\frac{d^2u}{dx} +3 = 0 $$ $$
 *
 $$ \displaystyle (Eq. 3.11.18)
 * }

From $$ Eq. 3.11.13 \quad $$, $$ u^h \quad $$ was defined as

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

$$ u^h = \sum_{i=0}^{n=2} d_ib_i = u^h(x) = d_0(1)   + d_1cos(x) +d_2sin(x) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.19)
 * }

Taking $$ 2^{nd} \quad $$ derivative with respect to "x"

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

$$ u^{h''} = d_00 + d_1(-cos(x)) + d_2(-sin(x)) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.20)
 * }

Therefore

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

$$ P(u^h) := 2(d_00 + d_1(-cos(x)) + d_2(-sin(x))) +3 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.21)
 * }

Afterwords, noting that

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

$$ \int_\omega w P(u^h) dx := 0 \Rightarrow w_i = b_i \Rightarrow \int_\omega b_i P(u^h) dx := 0 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.22)
 * }

which is also noted in $$ 10-4, Eq. (2) \quad $$

and substituting for $$ P(u^h) \quad $$ equation above, we get

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

 * style="width:95%" |
 * style="width:95%" |

$$ \int_\omega b_i (2(d_00 + d_1(-cos(x)) + d_2(-sin(x)) +3) dx := 0 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.11.23)
 * }

recalling that $$ \underline{b} =

\begin{bmatrix} 1  \\ cos(x) \\ sin(x) \end{bmatrix} $$

and taking any of the equations as the last needed equation, for example $$ b_1 \quad $$, we get

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

$$ \int_\omega
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} 1  \\ cos(x) \\ sin(x) \end{bmatrix}

[2(d_00 + d_1(-cos(x)) + d_2(-sin(x)) +3] dx := \int_\omega (-3

\begin{bmatrix} 1  \\ cos(x) \\ sin(x) \end{bmatrix} )dx \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.11.25)
 * }

Doing symbolic integration and evaluating interval on $$ \omega = ]0,1[ \quad $$ we get

The outputs are as follow:

$$ \int_\omega

\begin{bmatrix} 1  \\ cos(x) \\ sin(x) \end{bmatrix}

[1+cos(x)+sin(x)] dx = $$ 0 -1.683 -0.9194 0  -1.455 -0.7081 0 -0.7081 -0.5454

$$ = \int_\omega (-3

\begin{bmatrix} 1  \\ cos(x) \\ sin(x) \end{bmatrix} )dx \quad $$

-3.0 -2.524 -1.379

(4)

Assembling the matrix $$ \underline{K} \underline{d} = \underline{F} \quad $$ from 2 B.C. constrains and 1 more equation generated by weighting factor method:

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

$$ \begin{bmatrix} 1     & cos(1) & sin(1) \\ 0 &   -sin(0)  & cos(0) \\ 0 &  -1.683    & -0.9194     \end{bmatrix}
 * style="width:95%" |
 * style="width:95%" |

\begin{bmatrix} d_0     \\ d_1 \\ d_2 \end{bmatrix}

=

\begin{bmatrix} 0\\ -4 \\ -3   \end{bmatrix} $$ $$
 *
 $$ \displaystyle (Eq. 3.11.26)
 * }

For "K" matrix: 1.0 0.5403  0.8415  0      0      1.0  0   -1.683 -0.9194

For "F" matrix:

0 -4.0 -3.0

(5)

Solving for $$ \underline{d} = \underline{K}^{-1} F\quad $$ we get

d = 1.222 3.968 -4.0

(6)

Constructed $$ u^h \quad $$ matrix on matlab by using formula $$ u^h = \sum_{i=0}^{n=2} d_ib_i $$

3.968*cos(x) - 4.0*sin(x) + 1.222

Exact solution:

Given a general expression for 2nd Order ODE

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

$$ a_2(x)y'' + a_1(x)y + a_0(x)y = g(x) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.27)
 * }

and taking a more specific case ("Self adjoin") for this specific problem, we have

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

$$ [a_2(x)y']' + a_0(x)y = g(x) \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.17)
 * }

$$ g(x) = 0 \Rightarrow \quad $$ no forcing function (steady state)

$$ a_0(x)y = f(x) = 3 \quad $$ (given)

$$ a_2 = 2 \quad $$ (given)

Therefore $$ Eq. 3.11.17 \quad $$ becomes

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

 * style="width:95%" |
 * style="width:95%" |

$$ [2y']' + 3 = 0 \quad $$ $$
 *
 $$ \displaystyle (Eq. 3.11.28)
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ \frac{du}{dx}(\frac{du}{dx}2) + 3 = 0 $$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ \int \frac{du}{dx}2 = \int -3 dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle\!
 * }

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

$$ \frac{du}{dx}2 = -3x + c_1' \Rightarrow \frac{du}{dx} = -1.5x + c_1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.29)
 * }

at B.C of $$ -\frac{du}{dx}(0) = 4 \Rightarrow $$

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

 * style="width:95%" |
 * style="width:95%" |

$$ -4 = -1.5(0) + c_1 \Rightarrow c_1 = -4$$ $$
 *
 $$ \displaystyle\!
 * }

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

 * style="width:95%" |
 * style="width:95%" |

$$ \int {du} = \int (-1.5x - 4) dx $$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ u = \frac{-1.5}{2}x^2 - 4x + c_2 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.30)
 * }

at B.C of $$ u(1) = 0 \Rightarrow $$

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

 * style="width:95%" |
 * style="width:95%" |

$$ 0 = \frac{-1.5}{2}(1)^2 - 4(1) + c_2 \Rightarrow c_2 = -4.75$$ $$
 *
 $$ \displaystyle\!
 * }

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

$$ u(x) = -0.75x^2 - 4x - 4.75 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.11.31)
 * }

Circle view

So taken in account the exact solution formula, different plot were generated for different n values (n = 2, 4, 6)

Error came out to be:

x = 0.5

ex = n=2       n=4          n=6 -2.2405e-001 1.7223e-004 -1.9119e-004

Exact, approximation and error outputs

Approximation zoomed in view 251 equidistant point were taken for non-symbolic plot of n=2,4,6

--Egm6341.s11.team3.JA 02:32, 10 February 2011 (UTC)

= References =

= Contributing Members =

--Eml5526.s11.team5.srv 02:15, 13 February 2011 (UTC)

--Eml5526.s11.team5.smith 16:50, 15 February 2011 (UTC)

--Eml5526.s11.team5.savery 20:33, 15 February 2011 (UTC)

--Lokeshdahiya 14:19, 16 February 2011 (UTC)

--Eml5526.s11.team5.vijay 21:13, 16 February 2011 (UTC)

--Eml5526.s11.team5.mcdaniel 21:50, 2 February 2011 (UTC)

--Egm6341.s11.team3.JA 21:56, 16 February 2011 (UTC)

Eml5526.s11.team5.smith 02:09, 17 February 2011 (UTC)

title: User:Eml5526.s11.team5/sub3


text:
Colleges and Universities/University of Westminster

Wikiversity learning resources from courses at the University of Westminster
 * SPIR608 Political Simulations and Gaming

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


text: File:Plot of 'u'.png

title: File:Plot of 'u'.png


text: File:Zoom inEML5526s11.png

title: File:Zoom inEML5526s11.png


title: File:ErrorNMs11.png


text:
Template:Sec link/secure url/doc

This is the sec link/secure url sub-template.

Do not use this template directly, use sec link or sec link auto and their sister templates instead.

This template is called from the sec link family of templates.

Technical details
This sub-template returns a raw URL for the secure server. It expects these parameters:

For more documentation see sec link.

title: Template:Sec link/secure url/doc


text:
Template:Sec link/secure url

https://secure.wikimedia.org//wiki/

title: Template:Sec link/secure url