User:Eml5526.s11.team01/Hwk3

=Problem 1: Solving using WRF=

Solving a differential equation using Weight Residual Form. From lecture 13

Given
The partial differential equation is:


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

$$ \frac{\partial}{\partial x}\left(2\frac{\partial u}{\partial x}\right) + 3 = 0 $$
 * $$\displaystyle (Eq.1.1) $$
 * }
 * }

The boundary conditions are:

$$ u\left(1\right) = 0$$

$$\frac{du}{dx}\left(0\right) = -4 $$

A weighting funtion should be used:


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



$$ b_{i} = b_{i}=\left ( x+1 \right )^{i} $$
 * $$\displaystyle (Eq.1.2) $$
 * }
 * }

Find

 * 1.1 Find an approximate solution of the form $$U^{h}=\sum_{i=0}^{n}d_{i}b_{i}$$ for n=2
 * 1.2 Find two equations that enforce the boundary conditions
 * 1.3 Project the weight residues
 * 1.4 Display the equations in matrix form
 * 1.5 Solve for $$\displaystyle d $$
 * 1.6 Construct $$U^{h}$$ and plot $$U^{h}$$ and $$u$$
 * 1.7 Repeat 1.1 to 1.6 for n = 4 and n = 6

1.1 Find an approximate solution with n=2
When n = 2 we have:
 * {| style="width:100%" border="0"

$$ d_{i}=\begin{bmatrix} d0 & d1 & d2 \end{bmatrix} $$
 * $$\displaystyle (Eq.1.3) $$
 * }
 * }

and


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

$$ b_{i}=\begin{bmatrix} 1 & \left ( x+ 1\right) & \left ( x+ 1\right )^{2} \end{bmatrix} $$
 * $$\displaystyle (Eq.1.4) $$
 * }
 * }

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

$$ U^{h}=d0+d1\left ( x+1 \right )+d2\left ( x+1\right )^{2} $$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * $$\displaystyle (Eq. 1.5) $$
 * style= |
 * }

1.2 Find two equations that enforce the boundary conditions
The boundary conditions are:
 * {| style="width:100%" border="0"

$$ u\left(1\right)=0=d0+2d1+4d2$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * $$\displaystyle (Eq. 1.6) $$
 * style= |
 * }

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

$$ \frac{du}{dx}(0)=-4 = d1+2d2$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * $$\displaystyle (Eq. 1.7) $$
 * style= |
 * }

1.3 Project the weight residues
In this case changing the function u(x) by the approximated function $$U^{h}$$ in the PDE (Eq.1.1), we found that P(u) is:


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

$$ P\left ( U^{h} \right )=2\frac{\partial^2 U^{h}}{\partial x^2}+3=2d2+3 $$
 * $$\displaystyle (Eq.1.8) $$
 * }
 * }

The (Eq.1.8) can be written as:


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

$$ P\left ( U^{h} \right )= \begin{bmatrix}0 & 0 & 2\end{bmatrix}\begin{bmatrix} d0\\ d1\\

d2\end{bmatrix}+3 $$
 * $$\displaystyle (Eq.1.9) $$
 * }
 * }

Projecting the residue we have:


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

$$ \int_{a}^{b}b(x)P\left ( U^{h} \right )= 0 $$ After substitution of (Eq.1.4) and (Eq.1.9) in (Eq.1.10) the following equation is obtained
 * $$\displaystyle (Eq. 1.10) $$
 * }
 * }
 * {| style="width:100%" border="0"

$$ 2\int_{0}^{1}\begin{bmatrix} 1\\ \left (x+ 1\right ) \\ \left ( x+1\right )^{2} \end{bmatrix}\begin{bmatrix} 0 &0 & 2\end{bmatrix}dx\begin{bmatrix} d0\\ d1\\
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |

d2\end{bmatrix}=-3\int_{0}^{1} \begin{bmatrix} 1\\ \left ( x+1\right ) \\

\left ( x+1 \right )^{2}\end{bmatrix}dx $$


 * $$\displaystyle (Eq. 1.11) $$
 * style= |
 * }

1.4 Display the equations in matrix form
Performing the product and then the integration indicated in (Eq.1.11) we have:
 * {| style="width:100%" border="0"

$$ \begin{bmatrix} 0 &0 & 4 \\ 0 & 0& 6 \\ 0 & 0& 3\frac{28}{3}\end{bmatrix}\begin{bmatrix} d0\\ d1\\
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |

d2\end{bmatrix}=\begin{bmatrix} -3\\ -\frac{9}{2}\\

-21\end{bmatrix} $$


 * $$\displaystyle (Eq. 1.12) $$
 * style= |
 * }

The latter equation is clearly of the form $$\displaystyle Kd=F$$. We can observe that K is not symmetric.

1.5 Solve for $$ d $$
In (Eq.1.12) we have three equations; but, as we need to enforce the boundary conditions we take the only one of them and solve it together with (Eq.1.6) and (Eq.1.7). From these equations we now have the system:


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

$$ \begin{bmatrix} 1 &2   & 4 \\ 0 &1 &2 \\ 0 & 0& 4 \end{bmatrix}\begin{bmatrix} d0\\ d1\\

d2\end{bmatrix}=\begin{bmatrix} 0\\ -4\\

-3\end{bmatrix}$$
 * $$\displaystyle (Eq.1.13) $$
 * }
 * }

After solving the system in (Eq.1.13) by using a calculator, we obtain d:
 * {| style="width:100%" border="0"

$$\begin{bmatrix} d0\\ d1\\
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |

d2\end{bmatrix}=\begin{bmatrix} 8\\ -2.5\\ -0.75\end{bmatrix}$$


 * $$$$
 * style= |
 * }

1.6 Construct $$U^{h}$$ and plot $$U^{h}$$ and $$u$$
Replacing d in (Eq.1.5) we obtain the approximated solution:


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

$$ U^{h}=8-2.5\left ( x+1\right )-0.75\left ( x+1\right )^2 $$
 * $$\displaystyle (Eq.1.14) $$
 * }
 * }

On the other hand, the exact solution of the PDE (Eq.1.1) with the given boundary condition is:


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

$$ u\left ( x \right )=-\frac{3}{4}x^{2}-4x+\frac{19}{4} $$
 * $$\displaystyle (Eq.1.15) $$
 * }
 * }

We observe that the aproximated and the exact solution are the same; therefore, increasing the number of terms in b(x) will not improve the solution

In Fig.1.1 we show the exact and the approximated solution.



1.7 Repeat 1.1 to 1.6) for i = 4 and i = 6
Now we repeat the procedure using i = 4 and i=6 in Eq.1.4

As the procedure is the same shown before, we develop a routine using Matlab to perform the products, the integration and to solve the system of equations. The routine is show in the appendix. The error was also Zero for i = 4, and i = 6.

Appendix
 Matlab Code: 

=Problem 2: Spring System =

Given
Spring System with closed ends.

Find

 * 2.1 Number the Element and Nodes.
 * 2.2 Assemble the Global Stiffness and Force Matrix.
 * 2.3 Partition the system and solve for the Nodal Displacements.
 * 2.4 Compute the Reaction Forces.

2.1 Number the Element and Nodes.
The Spring Systems is shown in the figure.



Superscripts to each of the stiffness values denotes the element number. Furthermore, Node Numbers are also shown in the figure. As shown in the figure, there are 4 nodes.

2.2 Assemble the Global Stiffness and Force Matrix.
For Element Number 1,$$I=1 $$ and $$J=4$$,

$$K^{(1)}=k\begin{bmatrix} 1& -1\\ -1 & 1 \end{bmatrix}\Rightarrow \tilde{K}=\begin{bmatrix} k& 0 & 0 & -k\\ 0 & 0 & 0 & 0\\ 0& 0 & 0 & 0\\ -k& 0 & 0 & k \end{bmatrix}$$

For Element Number 2,$$I=4 $$ and $$J=3$$,

$$K^{(2)}=2k\begin{bmatrix} 1& -1\\ -1 & 1 \end{bmatrix}\Rightarrow \tilde{K}=\begin{bmatrix} 0& 0 & 0 & 0\\ 0 & 0 & 0 & 0\\ 0& 0 & 2k& -2k\\ 0& 0 & -2k & 2k \end{bmatrix}$$

For Element Number 3,$$I=3 $$ and $$J=2$$,

$$K^{(3)}=k\begin{bmatrix} 1& -1\\ -1 & 1 \end{bmatrix}\Rightarrow \tilde{K}=\begin{bmatrix} 0& 0 & 0 & 0\\ 0 & k & -k & 0\\ 0& -k & k& 0\\ 0& 0 & 0& 0 \end{bmatrix}$$

For Element Number 4,$$I=1 $$ and $$J=3$$,

$$K^{(4)}=3k\begin{bmatrix} 1& -1\\ -1 & 1 \end{bmatrix}\Rightarrow \tilde{K}=\begin{bmatrix} 3k& 0 & -3k & 0\\ 0 & 0& 0& 0\\ -3k& 0 & 3k& 0\\ 0& 0 & 0& 0 \end{bmatrix}$$

Assembled Global Systems Matrix will be given by,
 * {| style="width:100%" border="0"

$$K=\sum_{e=1}^{4}\tilde{K}^{e}=\begin{bmatrix} k+3k & 0 & -3k & -k\\ 0 & k & -k & 0\\ -3k & -k &6k &-2k \\ -k & 0 & -2k & 3k \end{bmatrix}$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style= |
 * }

Where Superscript $$'e'$$ denotes the Element Number.

The Displacement and Global Force Matrices for this system are given by,

2.3 Partition the system and solve for the Nodal Displacements.
Eventually, The Global System of Equations is Given by,

$$\begin{bmatrix} k+3k & 0 & -3k & -k\\ 0 & k & -k & 0\\ -3k & -k &6k &-2k \\ -k & 0 & -2k & 3k \end{bmatrix}\begin{bmatrix} 0\\ 0\\ u_{3}\\ u_{4} \end{bmatrix}=\begin{bmatrix} r_{1}\\ r_{2}\\ 50\\ 0 \end{bmatrix}$$ <p style="text-align:right;">$$\displaystyle (Eq.2.1) $$

As the first two displacements are prescribed, We can Partition the system after 2 Rows and 2 Columns,

And Can form the Matrix System of Equivalent to,

$$\begin{bmatrix} K_{E} & K_{EF}\\ K_{EF}^{T} & K_{F} \end{bmatrix}\begin{bmatrix} \bar{d}_{E}\\ d_{F} \end{bmatrix}=\begin{bmatrix} r_{E}\\ f_{F} \end{bmatrix}$$

Where,

$$K_{E}=\begin{bmatrix} 4k & 0\\ 0 & K \end{bmatrix}, K_{EF}=\begin{bmatrix} -3k & -k\\ -k& 0 \end{bmatrix},K_{EF} ^T=\begin{bmatrix} -3k & -k\\ -k & 0 \end{bmatrix}, K_{F}=\begin{bmatrix} 6k & -2k\\ -2k & 3K \end{bmatrix}, \bar{d}_{E}=\begin{bmatrix} 0\\ 0\\ \end{bmatrix}, d_{F}=\begin{bmatrix} u_{3}\\ u_{4} \end{bmatrix},r_{E}=\begin{bmatrix} r_{1}\\ r_{2} \end{bmatrix}, f_{F}=\begin{bmatrix} 50\\ 0 \end{bmatrix}$$

Solving the Following Reduced System of Equations,

$$K_{E}\bar{d}_{E}+K_{EF}d_{F}=r_{E}, K_{EF}^T\bar{d}_{E}+K_{F}d_{F}=f_{F}$$

Using $$\displaystyle (Eq.2.1) $$ we can form the Following Simultaneous Equations,

$$ 6k.u_{3}-2k.u_{4}=50$$ <p style="text-align:right;">$$\displaystyle (Eq.2.2) $$

$$ -2k.u_{3}-3k.u_{4}=0$$ <p style="text-align:right;">$$\displaystyle (Eq.2.3) $$

Solving These 2 Simultaneous Equations, We get,

2.4 Compute the Reaction Forces.
Now, Utilizing the First row of the math>\displaystyle (Eq.2.1) we can get

And Further, Utilizing the Second row of the math>\displaystyle (Eq.2.1) we can get,

=Problem 3: Truss Structure=

Given
The truss Structure Shown in the Figure. Nodes A & B are fixed. A Force of 10N Acts in the Positive x-Direction at node C. Young's Modulus is $$E=10^{11}Pa$$ and the Cross Sectional Area for all bars are $$ A= 2.10^{-2}$$

Find

 * 3.1 Number the Element and Nodes
 * 3.2 Assemble the Global Stiffness and Force Matrix
 * 3.3 Partition the system and solve for the Nodal Displacements
 * 3.4 Compute the Stresses & Reactions

3.1 Number the Element and Nodes
The Truss System Along with the Element Numbers and the Node Numbers Written on it, is Shown in the figure.



3.2 Assemble the Global Stiffness and Force Matrix
For Element Number 1,$$I=1 $$ and $$J=3$$, Also Element 1 is inclined at 90 deg to the Positive X-direction,$$\phi^{(1)} =90^{0}$$, $$cos(90)=0,sin(90)=1$$, $$l=1$$,

Using the Following formula for obtaining the Elemental Stiffness Matrix for all the elements,

$$K^{e}=\frac{AE}{l}\begin{bmatrix} (cos\phi ^{e})^{2} & cos\phi ^{e}.sin\phi ^{e} & -(cos\phi ^{e})^{2} & -cos\phi ^{e}.sin\phi ^{e} \\ cos\phi ^{e}.sin\phi ^{e} & (sin\phi ^{e})^{2} & -cos\phi ^{e}.sin\phi ^{e} &-(sin\phi ^{e})^{2} \\ -(cos\phi ^{e})^{2} & -cos\phi ^{e}.sin\phi ^{e} & (cos\phi ^{e})^{2} & cos\phi ^{e}.sin\phi ^{e}\\ -cos\phi ^{e}.sin\phi ^{e} & -(sin\phi ^{e})^{2} & cos\phi ^{e}.sin\phi ^{e} & (sin\phi ^{e})^{2} \end{bmatrix}$$ <p style="text-align:right;">$$\displaystyle (Eq.3.1) $$

Utilizing above equation for Element 1,

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

For Element Number 2,$$I=2$$ and $$J=4$$, Also Element 2 is inclined at 90 deg to the Positive X-direction,$$\phi^{(2)} =90^{0}$$, $$cos(90)=0,sin(90)=1$$, $$l=1$$, Therefore,

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

For Element Number 3,$$I=2$$ and $$J=3$$, Also Element 3 is inclined at 135 deg to the Positive X-direction,$$\phi^{(3)} =135^{0}$$, $$cos(135)=\frac{-1}{\sqrt{2}},sin(135)=\frac{1}{\sqrt{2}}$$, $$l=\sqrt{2}$$, Therefore,

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

For Element Number 4,$$I=3$$ and $$J=4$$, Also Element 4 is inclined at 90 deg to the Positive X-direction,$$\phi^{(2)} =0^{0}$$, $$cos(0)=1,sin(0)=0$$, $$l=1$$, Therefore,

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

Making the Assembly of all these 4 Elemental Matrices, we get the following Global Stiffness Mtria for the given Truss System,

<p style="text-align:right;">$$\displaystyle (Eq.3.2) $$ Now, Only force is acting in positive X-Direction at C, and there will be Reactions at A & B.

Therefore, The Force Matrix Can be Assembled as,

<p style="text-align:right;">$$\displaystyle (Eq.3.3) $$

3.3 Partition the system and solve for the Nodal Displacements
From $$\displaystyle (Eq.3.2) $$ & $$\displaystyle (Eq.3.3)  $$ The Global System os Equation will be,

$$AE \begin{bmatrix} 0 & 0 &0 & 0 & 0 & 0 &  0& 0\\ 0 & 1 & 0 & 0 & 0 & -1 & 0 &0 \\ 0 & 0 & \frac{1}{2\sqrt{2}} & \frac{-1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} &0  &0 \\ 0 & 0 & -\frac{1}{2\sqrt{2}} &1+\frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & 0 & -1\\ 0 & 0 & -\frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} & 1+\frac{1}{2\sqrt{2}} &-\frac{1}{2\sqrt{2}} & -1 &0 \\ 0 &-1 &\frac{1}{2\sqrt{2}}  & -\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} &1+\frac{1}{2\sqrt{2}}  & 0 & 0\\ 0 & 0 & 0 & 0 & -1 & 0 & 1 & 0\\ 0 & 0 & 0 & -1 & 0 & 0 & 0 & 1 \end{bmatrix}.\begin{bmatrix} 0\\ 0\\ 0\\ 0\\ u_{3x}\\ u_{3y}\\ u_{4x}\\ u_{4y} \end{bmatrix}=\begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y}\\ 0\\ 0\\ 10\\ 0 \end{bmatrix}$$

As the Displacements are prescribed untill second node, Partitioning the systems after second row and column and further comparing with the following equation,

$$\begin{bmatrix} K_{E} & K_{EF}\\ K_{EF}^{T} & K_{F} \end{bmatrix}\begin{bmatrix} \bar{d}_{E}\\ d_{F} \end{bmatrix}=\begin{bmatrix} r_{E}\\ f_{F} \end{bmatrix}$$

we get, $$\bar{d_{E}}=\begin{bmatrix} 0\\ 0\\ 0\\ 0 \end{bmatrix},\ d_{F}=\begin{bmatrix} u_{3x}\\ u_{3y}\\ u_{4x}\\ u_{4y} \end{bmatrix},\ r_{E}=\begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y} \end{bmatrix},\ f_{F}=\begin{bmatrix} 0\\ 0\\ 10\\ 0 \end{bmatrix},\ K_{F}=\begin{bmatrix} 1+\frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} & -1 &0 \\ -\frac{1}{2\sqrt{2}} &1+\frac{1}{2\sqrt{2}} &0  & 0\\ -1 & 0 & 1 & 0\\ 0 & 0 & 0 &1 \end{bmatrix},\ K_{EF}=\begin{bmatrix} 0 & 0 & 0 & 0\\ 0 &-1 &0  & 0\\ -\frac{1}{2\sqrt{2}} & \frac{1}{2\sqrt{2}} &  0& 0\\ \frac{1}{2\sqrt{2}} & -\frac{1}{2\sqrt{2}} &0 &-1 \end{bmatrix}$$

Eventually, we can obtain,

$$d_{F}=K_{F}^{-1}(f_{F}-K_{EF}^T.\bar{d}_{E})$$

But $$\bar{d}_{E}$$ is a Zero Matrix, which will yield, $$d_{F}=K_{F}^{-1}.f_{F}$$ Solving this system in MATLAB, Using the following code, we get,

 Matlab Code: 

<p style="text-align:right;">$$\displaystyle (Eq.3.4) $$

These are the Unknown Nodal Displacements.

3.4 Compute the Stresses & Reactions
Now the Reactions can be calculated as,

$$r_{E}=K_{E}.\bar{d}_{E}+K_{EF}.d_{F}$$

But again $$\bar{d}_{E}$$ is a Zero Matrix, Which will give, $$r_{E}=K_{EF}.d_{F}$$.

Using the Following MATLAB Code to obtain the Reaction Matrix,

 Matlab Code: 

<p style="text-align:right;">$$\displaystyle (Eq.3.5) $$

These are the Unknown Reactions. Now We can calculate the Stress in each element using the following formula,

$$\sigma ^e=\frac{E^e}{l^e}\begin{bmatrix} -cos\phi ^{e} &-sin\phi ^{e} & cos\phi ^{e} & sin\phi ^{e} \end{bmatrix}.d^e$$ <p style="text-align:right;">$$\displaystyle (Eq.3.6) $$

Now we can utilize the above equation to calculate the Stresses for all the elements,

For element 1, $$\phi =90^0, l=1$$ &

$$d^1=\begin{bmatrix} u_{1x}\\ u_{1y}\\ u_{3x}\\ u_{3y} \end{bmatrix}=\frac{1}{2\times 10^{9}}\begin{bmatrix} 0\\ 0\\ 38.2843\\ 10 \end{bmatrix}$$

Using the Following MATLAB code,  Matlab Code: 

Therefore,

In the similar fashion, For Element 2, Utilizing $$\displaystyle (Eq.3.6) $$ we get,

For Element 3, Utilizing $$\displaystyle (Eq.3.6) $$ we get,

Finally, For Element 4, Utilizing $$\displaystyle (Eq.3.6) $$ we get,

=Problem 4 : Show that $$\displaystyle \alpha$$ does not equal $$\displaystyle   \beta$$ =

Given
Given


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

$$\alpha =b_{i}\left[ a_{2}^{\prime }b_{j}^{\prime }+a_{2}b_{j}^{\prime \prime } \right]$$ And
 * <p style="text-align:right;">$$\displaystyle (Eq. 3.1) $$
 * }
 * }


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



$$\beta =b_{j}\left[ a_{2}^{\prime }b_{i}^{\prime }+a_{2}b_{i}^{\prime \prime } \right]$$
 * <p style="text-align:right;">$$\displaystyle (Eq. 3.2) $$
 * }
 * }

Find
Show that $$\alpha \ne \beta $$

Solution
Let


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

$$b_{i}\left( x \right)=\cos \left( ix \right)$$ And
 * <p style="text-align:right;">$$\displaystyle (Eq. 3.3) $$
 * }
 * }
 * {| style="width:100%" border="0"

$$b_{j}\left( x \right)=\cos \left( jx \right)$$ Assume that $$i=1$$ and $$j=2$$. Therefore, equations 3.3 and 3.4 become:
 * <p style="text-align:right;">$$\displaystyle (Eq. 3.4) $$
 * }
 * }
 * {| style="width:100%" border="0"

$$b_{i}\left( x \right)=\cos \left( x \right)$$
 * <p style="text-align:right;">$$$$
 * }
 * }


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

$$b_{j}\left( x \right)=\cos \left( 2x \right)$$ Differentiating:
 * <p style="text-align:right;">$$$$
 * }
 * }
 * {| style="width:100%" border="0"

$$b_{i}^{\prime }\left( x \right)=-\sin \left( x \right)$$
 * <p style="text-align:right;">$$$$
 * }
 * }


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

$$b_{j}^{\prime }\left( x \right)=-2\sin \left( 2x \right)$$
 * <p style="text-align:right;">$$$$
 * }
 * }


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

$$b_{i}^{\prime \prime }\left( x \right)=-\cos \left( x \right)$$
 * <p style="text-align:right;">$$$$
 * }
 * }


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

$$b_{j}^{\prime \prime }\left( x \right)=-4\cos \left( 2x \right)$$ Plugging in values and re-arranging
 * <p style="text-align:right;">$$$$
 * }
 * }
 * {| style="width:100%" border="0"

$$\alpha =-2a_{2}^{\prime }\sin \left( 2x \right)\cos \left( x \right)-4a_{2}^{\prime \prime }\cos \left( 2x \right)\cos \left( x \right)$$
 * <p style="text-align:right;">$$\displaystyle (Eq. 3.5) $$
 * }
 * }


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

$$\beta =-a_{2}^{\prime }\sin \left( x \right)\cos \left( 2x \right)-a_{2}^{\prime \prime }\cos \left( x \right)\cos \left( 2x \right)$$ By inspection of terms $$\displaystyle a_{2}^{\prime }$$ and $$\displaystyle a_{2}^{\prime \prime }$$ it can be shown that
 * <p style="text-align:right;">$$\displaystyle (Eq. 3.6) $$
 * }
 * }
 * {| style="width:100%" border="0"

$$\alpha \ne \beta $$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * <p style="text-align:right;">$$\displaystyle ( Eq. 3.7) $$
 * style= |
 * }

=Problem 5: Equivalent Stiffness of a Spring =

Given
A spring is given in Figure 5.1 as follows:

Find
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 Figure 5.1 is can be described by


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

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


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.1) $$
 * }
 * }

Solution
The left end is constrained, i.e., prescribed displacement is zero at the left end.

Assume there is a force F acting on the right end, and the nodes are numbered in (Fig.5.2).



The Element stiffness matrices are the following:


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

$$ \mathbf{K}^{\left ( 1\right )}  =\begin{bmatrix} k^{\left ( 1 \right )} & -k^{\left ( 1 \right )} \\ -k^{\left ( 1 \right )} & k^{\left ( 1 \right )} \end{bmatrix},

\mathbf{K}^{\left ( 2\right )}  =\begin{bmatrix} k^{\left ( 2 \right )} & -k^{\left ( 2 \right )} \\ -k^{\left ( 2 \right )} & k^{\left ( 2 \right )} \end{bmatrix},

\mathbf{K}^{\left ( 3\right )}  =\begin{bmatrix} k^{\left ( 3 \right )} & -k^{\left ( 3 \right )} \\ -k^{\left ( 3 \right )} & k^{\left ( 3 \right )} \end{bmatrix}$$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.2)$$
 * }
 * }

By direct assembly, the global stiffness matrix is:


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

$$ \mathbf{K}=\begin{bmatrix} k^{\left ( 1 \right )} & -k^{\left ( 1 \right )} & 0 & 0\\ -k^{\left ( 1 \right )} & k^{\left ( 1 \right )} + k^{\left ( 2 \right )} & -k^{\left ( 2 \right )} & 0\\ 0 & -k^{\left ( 2 \right )}  & k^{\left ( 2 \right )} + k^{\left ( 3 \right )}  & -k^{\left ( 3 \right )}   \\ 0 & 0 & -k^{\left ( 3 \right )}  & k^{\left ( 3 \right )} \end{bmatrix} $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.3) $$
 * }
 * }

The displacements and force matrices for the system are:


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

$$ \mathbf{d}=\begin{bmatrix} 0\\ u_{2}\\ u_{3}\\ u_{4} \end{bmatrix},

\mathbf{f}=\begin{bmatrix} 0\\ 0\\ 0\\ F \end{bmatrix},

\mathbf{r}=\begin{bmatrix} r_{1}\\ 0\\ 0\\ 0\end{bmatrix} $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.4) $$
 * }
 * }

The global system of equations is given by:


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

$$ \begin{bmatrix} k^{\left ( 1 \right )} & -k^{\left ( 1 \right )} & 0 & 0\\ -k^{\left ( 1 \right )} & k^{\left ( 1 \right )} + k^{\left ( 2 \right )} & -k^{\left ( 2 \right )} & 0\\ 0 & -k^{\left ( 2 \right )}  & k^{\left ( 2 \right )} + k^{\left ( 3 \right )}  & -k^{\left ( 3 \right )}   \\ 0 & 0 & -k^{\left ( 3 \right )}  & k^{\left ( 3 \right )} \end{bmatrix}

\begin{bmatrix} 0\\ u_{2}\\ u_{3}\\ u_{4} \end{bmatrix}=

\begin{bmatrix} r_{1}\\ 0\\ 0\\ F \end{bmatrix} $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.5) $$
 * }
 * }

Solve the matrix equation above; we can have the following four equations


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

$$ \left\{\begin{matrix} -k^{\left ( 1 \right )}u_{2}=r_{1} \\ u_{2}\left ( k^{\left (1 \right )} + k^{\left (2 \right)} \right )- u_{3}k^{\left (2 \right)}=0 \\ -u_{2}k^{\left (2 \right )} + u_{3}\left ( k^{\left (2 \right )} + k^{\left (3 \right)} \right ) - u_{4}k^{\left (3 \right)}=0  \\ -u_{3}k^{\left (3 \right)}+u_{4}k^{\left (3 \right)}=F \end{matrix}\right. $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.6) $$
 * }
 * }

Because we have the following relationship:


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

$$ k=\frac{AE}{l} $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.7) $$
 * }
 * }

Then we can get the following results:
 * {| style="width:100%" border="0"

$$ k^{\left ( 1 \right )}=\frac{10atE}{l}, k^{\left ( 2 \right )}=\frac{5btE}{l}, k^{\left ( 3 \right )}=\frac{10atE}{l} $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.8) $$
 * }
 * }

Substitute (Eq.5.8) into (Eq.5.6), we can have the following result:


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

$$ u_{4}=\frac{Fl\left ( a+b \right )}{5abtE} $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.9) $$
 * }
 * }

where $$u_{4}$$ is the displacement of node 4.

And then,


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

$$ \frac{F}{u_{4}}=\frac{5Etab}{\left ( a+b \right )l} $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 5.10) $$
 * }
 * }

which yields:


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

$$ k=\frac{5Etab}{\left ( a+b \right )l} $$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * <p style="text-align:right;">$$\displaystyle ( Eq. 5.11) $$
 * style= |
 * }

=Problem 6: Three Bar Structure =

Given
Given a three bar structure subjected to the prescribed load at point C equal to $$ 10^3 ~N $$ as shown in Figure 2.19 of the book. The Young's modulus is $$ E = 10^{11} ~Pa $$, the cross-sectional area of the bar BC is $$ 2 \times 10^{-2} ~m^2 $$ and that of BD and BF is $$ 10^{-2} ~m^2 $$. Note that point D is free to move in the x-direction. Coordinates of joints are given in meters.



Find

 * 6.1 Construct the global stiffness matrix and load matrix.


 * 6.2 Partition the matrices and solve for the unknown displacements at Point B and displacement in the x-direction at point D.


 * 6.3 Find stresses in the three bars.


 * 6.4 Find the reactions at nodes C,D,F.

6.1 Construct the global stiffness matrix and load matrix
For Element 3, $$ \phi = 45^o,I = 3,J =4 $$

$$ K^{\left(3\right)} = \frac{10^{-2} \times 10^{11}}{\sqrt{2}} \begin{bmatrix} \frac{1}{2} &\frac{1}{2} &-\frac{1}{2} &-\frac{1}{2}\\ \frac{1}{2} &\frac{1}{2} &-\frac{1}{2} &-\frac{1}{2}\\ -\frac{1}{2} &-\frac{1}{2} &\frac{1}{2} &\frac{1}{2}\\ -\frac{1}{2} &-\frac{1}{2} &\frac{1}{2} &\frac{1}{2} \end{bmatrix} $$

For Element 2, $$ \phi = 90^o,I = 2,J = 4 $$

$$ K^{\left(2\right)} = \frac{2 \times 10^{-2} \times 10^{11}}{1} \begin{bmatrix} 0 &0  &0  &0\\  0  &1  &0 &-1\\  0  &0  &0 &0 \\  0  &-1  &0 &1 \end{bmatrix} $$

For Element 1, $$ \phi = 135^o,I = 1,J =4 $$

$$ K^{\left(1\right)} = \frac{10^{-2} \times 10^{11}}{\sqrt{2}} \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} $$

Global Stiffness Matrix

$$ K = 10^{9} \begin{bmatrix} \frac{1}{2 \sqrt{2}} &-\frac{1}{2 \sqrt{2}} &0  &0  &0  &0  &-\frac{1}{2 \sqrt{2}}  &\frac{1}{2 \sqrt{2}} \\ -\frac{1}{2 \sqrt{2}} &\frac{1}{2 \sqrt{2}} &0  &0  &0  &0  &\frac{1}{2 \sqrt{2}}  &-\frac{1}{2 \sqrt{2}} \\ 0 &0 &0  &0  &0  &0  &0  &0 \\ 0 &0  &0  &2  &0  &0  &0  &-2 \\ 0 &0  &0  &0  &\frac{1}{2 \sqrt{2}}  &\frac{1}{2 \sqrt{2}}  &-\frac{1}{2 \sqrt{2}}  &-\frac{1}{2 \sqrt{2}} \\ 0 &0 &0  &0  &\frac{1}{2 \sqrt{2}}  &\frac{1}{2 \sqrt{2}}  &-\frac{1}{2 \sqrt{2}}  &-\frac{1}{2 \sqrt{2}} \\ -\frac{1}{2 \sqrt{2}} &\frac{1}{2 \sqrt{2}} &0  &0  &-\frac{1}{2 \sqrt{2}}  &-\frac{1}{2 \sqrt{2}}  &\frac{1}{\sqrt{2}}  &0 \\ \frac{1}{2 \sqrt{2}} &-\frac{1}{2 \sqrt{2}} &0  &-2  &-\frac{1}{2 \sqrt{2}}  &-\frac{1}{2 \sqrt{2}}  &0  &\frac{1}{\sqrt{2}}+2 \end{bmatrix} $$

$$ d = \begin{bmatrix} 0\\ 0\\ 0\\ 0\\ u_{3x}\\ 0\\ u_{4x}\\ u_{4y} \end{bmatrix} f = \begin{bmatrix} 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 10^3\\ 0 \end{bmatrix} r = \begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y}\\ 0\\ r_{3y}\\ 0\\ 0 \end{bmatrix} $$

Therefore

$$ load~matrix = \begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y}\\ 0\\ r_{3y}\\ 10^3\\ 0 \end{bmatrix} $$

6.2 Partition the matrices and solve
Global system of equations:

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

Partitions

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

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

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

f_{F} = \begin{bmatrix} 0\\ r_{3y}\\ 10^3\\ 0 \end{bmatrix} $$

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

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

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

K_{E} = \begin{bmatrix} \frac{1}{2 \sqrt{2}} &-\frac{1}{2 \sqrt{2}} &0  &0 \\ -\frac{1}{2 \sqrt{2}} &\frac{1}{2 \sqrt{2}} &0  &0 \\ 0 &0 &0  &0 \\ 0 &0  &0  &2 \end{bmatrix} $$

Unknown displacement at B is solved by

$$ K_{EF}^{T}\bar{d}_{E} + K_{F}d_{F} = f_{F} $$

Subtracting the first term from both sides of the above equation and premultiply by $$ K_F^{-1} $$, we obtain

$$ d_F = K_F^{-1} \left( f_F - K_{EF}^T \bar{d}_E \right) $$

Solving for the displacements we get: $$ d_{F}= 10^{-6} \begin{bmatrix} 1+2\sqrt{2}\\ 0\\ \frac{1}{2}+2\sqrt{2}\\ \frac{1}{2} \end{bmatrix} ~m $$

6.3 Find stresses in the three bars
Element 3 Stress

$$ \sigma = 0 ~Pa $$

Element 2 Stress

$$ \sigma = \frac{E}{l} \begin{bmatrix} 0 &-1 &0  &1 \end{bmatrix} 10^{-6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 0\\ 0 \end{bmatrix} = -5 \times 10^4 ~Pa $$

Element 1 Stress

$$ \sigma = \frac{E}{l} \begin{bmatrix} \frac{\sqrt{2}}{2} &-\frac{\sqrt{2}}{2} &\frac{\sqrt{2}}{2}  &\frac{\sqrt{2}}{2} \end{bmatrix} 10^{-6} \begin{bmatrix} \frac{1}{2} + 2\sqrt{2}\\ \frac{1}{2}\\ 0\\ 0 \end{bmatrix} = \sqrt{2} \times 10^5 ~Pa $$

6.4 Find the reactions at nodes C,D,F.
Reactions can be solved using this equation and the global matrix

$$ r_E = K_E \bar{d}_E + K_{EF} d_F $$

$$ load~matrix = \begin{bmatrix} r_{1x}\\ r_{1y}\\ r_{2x}\\ r_{2y}\\ 0\\ r_{3y}\\ 10^3\\ 0 \end{bmatrix}

= \begin{bmatrix} -10^3\\ 10^3\\ 0\\ -10^3\\ 0\\ 0\\ 10^3\\ 0 \end{bmatrix} ~ N $$

=Problem 7: Solving using WRF= Solving a differential equation using Weight Residual Form. From lecture 17

Given
The partial differential equation is:


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

$$ \frac{\partial}{\partial x}\left(2\frac{\partial u}{\partial x}\right) + 3 = 0 $$
 * <p style="text-align:right;">$$\displaystyle (Eq. 7.1) $$
 * }
 * }

The boundary conditions are:

$$ u\left(1\right) = 0$$

$$\frac{du}{dx}\left(0\right) = -4 $$

As weighting function should be used:


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



$$ b_{i} = b_{i}=\left ( x \right )^{i} $$
 * <p style="text-align:right;">$$\displaystyle (Eq. 7.2) $$
 * }
 * }

Find

 * 7.1 Find an approximate solution of the form $$U^{h}=\sum_{i=0}^{n}d_{i}b_{i}$$ for n=2
 * 7.2 Find two equations that enforce the boundary conditions
 * 7.3 Project the weight residues
 * 7.4 Display the equations in matrix form
 * 7.5 Solve for $$\displaystyle d $$
 * 7.6 Construct $$U^{h}$$ and plot $$U^{h}$$ and $$u$$
 * 7.7 Repeat 7.1 to 7.6 for n = 4 and n = 6

7.1 Find an approximate solution with n=2
When n = 2 we have:
 * {| style="width:100%" border="0"

$$d_{i}=\left[ \begin{matrix} d_{0} & d_{1} & d_{2} \\ \end{matrix} \right]$$
 * <p style="text-align:right;">$$\displaystyle (Eq.7.3) $$
 * }
 * }

and


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

$$b_{i}=\left[ \begin{matrix} 1 & \left( x \right) & \left( x \right)^{2} \\ \end{matrix} \right]$$
 * <p style="text-align:right;">$$\displaystyle (Eq.7.4) $$
 * }
 * }

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

$$U^{h}=d_{0}+d_{1}\left( x \right)+d_{2}\left( x \right)^{2}$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * <p style="text-align:right;">$$\displaystyle (Eq. 1.5) $$
 * style= |
 * }

7.2 Find two equations that enforce the boundary conditions
The boundary conditions are:
 * {| style="width:100%" border="0"

$$u\left( 1 \right)=0=d_{0}+d_{1}+d_{2}$$ Differentiating
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * <p style="text-align:right;">$$\displaystyle (Eq. 7.6) $$
 * style= |
 * }
 * {| style="width:100%" border="0"

$$\frac{du}{dx}=0+d_{1}+2d_{2}x$$
 * <p style="text-align:right;">$$$$
 * }
 * }

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


 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |

$$\frac{du}{dx}(0)=-4=d_{1}$$


 * <p style="text-align:right;">$$\displaystyle (Eq. 7.7) $$
 * style= |
 * }

7.3 Project the weight residues
In this case changing the function $$u\left( x \right)$$ by the approximated function $$U^{h}$$ in the PDE (Eq.7.1), we found that $$P\left( u \right)$$ is:


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

$$P\left( U^{h} \right)=2\frac{\partial ^{2}U^{h}}{\partial x^{2}}+3=4d_{2}+3$$
 * <p style="text-align:right;">$$\displaystyle (Eq.7.8) $$
 * }
 * }

The (Eq.7.8) can be written as:


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



$$P\left( U^{h} \right)=\left[ \begin{matrix} 0 & 0 & 4 \\ \end{matrix} \right]\left[ \begin{matrix} d_{0} \\ d_{1} \\ d_{2} \\ \end{matrix} \right]+3$$
 * <p style="text-align:right;">$$\displaystyle (Eq.7.9) $$
 * }
 * }

Projecting the residue we have:


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

$$ \int_{a}^{b}b(x)P\left ( U^{h} \right )= 0 $$ After substitution of equations 7.4 and 7.9) in 7.10, the following equation is obtained
 * <p style="text-align:right;">$$\displaystyle (Eq. 7.10) $$
 * }
 * }
 * {| style="width:100%" border="0"

$$\int_{0}^{1}{\left[ \begin{matrix} 1 \\   x  \\ x^{2} \\ \end{matrix} \right]}\left[ \begin{matrix} 0 & 0 & 4 \\ \end{matrix} \right]dx\left[ \begin{matrix} d_{0} \\ d_{1} \\ d_{2} \\ \end{matrix} \right]=-3\int_{0}^{1}{\left[ \begin{matrix} 1 \\   x  \\ x^{2} \\ \end{matrix} \right]}dx$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * <p style="text-align:right;">$$\displaystyle (Eq. 7.11) $$
 * style= |
 * }

7.4 Display the equations in matrix form
Performing the product and then the integration indicated in equation 7.11, we have
 * {| style="width:100%" border="0"


 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |

$$\left[ \begin{matrix} 0 & 0 & 4 \\   0 & 0 & 2  \\   0 & 0 & \frac{4}{3}  \\ \end{matrix} \right]\left[ \begin{matrix} d_{0} \\ d_{1} \\ d_{2} \\ \end{matrix} \right]=\left[ \begin{matrix} -3 \\   -\frac{3}{2}  \\ -1 \\ \end{matrix} \right]$$
 * <p style="text-align:right;">$$\displaystyle (Eq. 7.12) $$
 * style= |
 * }

The latter equation is clearly of the form $$\displaystyle Kd=F$$. We can observe that matrix K is not symmetric.

7.5 Solve for $$\displaystyle d $$
In equation 7.12 we have three equations; but, as we need to enforce the boundary conditions we take the only one of them and solve it together with equations 7.6 and 7.7. From these equations we now have the system:


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



$$\left[ \begin{matrix} 1 & 1 & 1 \\   0 & 1 & 2  \\   0 & 0 & 2  \\ \end{matrix} \right]\left[ \begin{matrix} d_{0} \\ d_{1} \\ d_{2} \\ \end{matrix} \right]=\left[ \begin{matrix} 0 \\   -4  \\   -3  \\ \end{matrix} \right]$$
 * <p style="text-align:right;">$$\displaystyle (Eq.7.13) $$
 * }
 * }

The latter is a linear system and can be solved by using any calculator. After solving the system in (Eq.7.13) we obtain d:
 * {| style="width:100%" border="0"


 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |

$$\left[ \begin{matrix} d_{0} \\ d_{1} \\ d_{2} \\ \end{matrix} \right]=\left[ \begin{matrix} \frac{5}{2} \\ -1 \\   -\frac{4}{3}  \\ \end{matrix} \right]$$
 * <p style="text-align:right;">$$$$
 * style= |
 * }

7.6 Construct $$U^{h}$$ and plot $$U^{h}$$ and $$u$$
Replacing $$d_{i}$$ in equation 7.5 we obtain the approximated solution:


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



$$U^{h}=\frac{5}{2}-x-\frac{4}{3}x^{2}$$
 * <p style="text-align:right;">$$\displaystyle (Eq.7.14) $$
 * }
 * }

On the other hand, the exact solution of the PDE (Eq.7.1) with the given boundary condition is:


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

$$ u\left ( x \right )=-\frac{3}{4}x^{2}-4x+\frac{19}{4} $$
 * <p style="text-align:right;">$$\displaystyle (Eq.7.15) $$
 * }
 * }

We observe that the approximated and the exact solution are not the same. Perhaps they would be closer if we increased the number of terms. Figure 7.1 we show the exact and the approximated solution.



7.7 Repeat 7.1 to 7.6) for i = 4 and i = 6
Now we repeat the procedure using n = 4 and n=6 in Eq.7.4

As the procedure is the same shown before, we develop a routine using Matlab to perform the products, the integration and to solve the system of equations. The routine is shown below in the appendix. Figures 7.2 and 7.3 show the results of n=4 and n=6, respectively. Figure 7.4 shows the error versus the number of terms.



Appendix
 Matlab Code: 

=Problem 8: Weak Form Proof=

Given
The Strong form is:
 * {| style="width:100%" border="0"

$$\displaystyle \frac{\mathrm{d} }{\mathrm{d} x} \left ( AE\frac{\mathrm{d} u}{\mathrm{d} x} \right ) + 2x=0$$ on      1<x<3,
 * <p style="text-align:right;">$$\displaystyle (Eq. 8.1a)  $$
 * }
 * }


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

$$ \displaystyle \sigma (1)=\left ( E\frac{\mathrm{d} u}{\mathrm{d} x} \right )_{x=1}=0.1, $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 8.1b)  $$
 * }
 * }


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

$$ \displaystyle u(3)=0.001 $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 8.1c)  $$
 * }
 * }

Find
Show that the weak form is:
 * {| style="width:100%" border="0"

$$ \int_{1}^{3}\frac{\mathrm{d} w}{\mathrm{d} x}AE\frac{\mathrm{d} u}{\mathrm{d} x}dx = -0.1\left ( wA \right )_{x=1}+\int_{1}^{3}2xwdx,\forall w$$ with $$w(3)=0$$


 * <p style="text-align:right;">$$\displaystyle (Eq. 8.2) $$
 * }
 * }

Solution
Equation (8.1b) is a condition on the derivative of u(x), so it is a natural boundary conditon; (8.1c) is a condition on u(x), so it is an essential boundary condition. Therefore, as the weight function must vanish on the essential boundaries, we consider all smooth weight functions w(x) such that w(3)=0. The trial solutions must satisfy the essential boundary condition u(3) = 0.001.

We start by multiplying the governing equation and the natural boundary condition over the domains where they hold by an arbitrary weight function:
 * {| style="width:100%" border="0"

$$ \int_{1}^{3} w [\frac{\mathrm{d} }{\mathrm{d} x}(AE\frac{\mathrm{d} u}{\mathrm{d} x})+2x]dx = 0, \forall w(x), $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 8.3) $$
 * }
 * }


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

$$ [wA(E\frac{\mathrm{d} u}{\mathrm{d} x}-0.1)]_{x=1} = 0, \forall w(1). $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 8.4) $$
 * }
 * }

Next we integrate (Eq.8.3) by parts, then we get the following result:


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

$$ \int_{1}^{3} [w (\frac{\mathrm{d} }{\mathrm{d} x}AE\frac{\mathrm{d} u}{\mathrm{d} x})]dx = \left (wAE\frac{\mathrm{d} u}{\mathrm{d} x} \right )\mid _{x=1}^{x=3} - \int_{1}^{3}  \frac{\mathrm{d} w }{\mathrm{d} x}AE\frac{\mathrm{d} u}{\mathrm{d} x}dx $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 8.5) $$
 * }
 * }

We have constructed the weight functions so that w(3)=0; therefore, the first term on the RHS of (Eq.8.5) vanishes at x=3. Substituting (Eq.8.5) into (Eq.8.3) gives:


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

$$ -\left (wAE\frac{\mathrm{d} u}{\mathrm{d} x} \right )\mid _{x=1} - \int_{1}^{3}  \frac{\mathrm{d} w }{\mathrm{d} x}AE\frac{\mathrm{d} u}{\mathrm{d} x}dx + \int_{1}^{3} 2xwdx=0, \forall w(x) $$ with w(3)=0. Substituting (Eq.8.4) into the first term of (Eq.8.6) gives:
 * <p style="text-align:right;">$$\displaystyle (Eq. 8.6) $$
 * }
 * }


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

$$ \int_{1}^{3} \frac{\mathrm{d} w }{\mathrm{d} x}AE\frac{\mathrm{d} u}{\mathrm{d} x}dx = -0.1(wA)_{x=1} + \int_{1}^{3} 2xwdx, \forall w(x) $$ with w(3)=0.
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * <p style="text-align:right;">$$\displaystyle ( Eq. 8.7) $$
 * style= |
 * }

=Problem 9: A solution to the weak form in Problem 3.1=

Consider a trial (candidate) solution of the form $$ u\left ( x \right ) = \alpha_{0} + \alpha_{1}\left ( x-3 \right ) $$ and a weight function of the same form. Obtain a solution to the weak form in Problem 3.1. Check the equilibrium equation in the strong form in Problem 3.1; is it satisfied?

Check the natural boundary condition; is it satisfied?

Solution
$$ u\left ( x \right ) = \alpha_{0} + \alpha_{1}\left ( x-3 \right ) $$

$$ w\left ( x \right ) = \beta_{0} + \beta_{1}\left ( x-3 \right ) $$

Using boundary conditions from Problem 3.9

$$ w \left(3 \right) = 0 \therefore  \beta_{0} = 0 $$

$$ u \left(3 \right) = \alpha_{0} = 0.001 $$

Only one unknown parameter and one arbitrary parameter remain

$$ u\left ( x \right ) = 0.001 + \alpha_{1}\left ( x-3 \right ) \frac{du}{dx} = \alpha_{1} $$

$$ w\left( x \right) = \beta_{1} \left(x - 3 \right) \frac{dw}{dx} = \beta_{1} $$

Now substitute into Problem 3.1 weak form

$$ \int_{1}^{3} \frac{\mathrm{d} w }{\mathrm{d} x}AE\frac{\mathrm{d} u}{\mathrm{d} x}dx =  -0.1(wA)_{x=1} + \int_{1}^{3} 2xwdx$$

$$ \int_{1}^{3} \beta_{1}\alpha_{1}EAdx + 0.1(\beta_{1} \left(x-3 \right)A)_{x=1} - \int_{1}^{3} 2\beta_{1}x \left(x-3 \right)dx = 0 $$

Integrating and factoring out $$ \beta_1 $$ we get:

$$ \beta_{1} \left( 2 \alpha_{1} EA - 0.2A + \frac{20}{3} \right) = 0 $$

The $$ \beta_1 $$ goes to zero and $$ \alpha_1 $$ can be solved:

$$ \alpha_{1} = \frac{0.1A - \frac{10}{3}}{EA} $$

Therefore the solution to the weak form is

$$ u \left( x \right) = 0.001 + \frac{0.1A - \frac{10}{3}}{EA} \left( x-3 \right) $$

Checking for equilibrium in the strong form, values are substituted into

$$ \frac{d}{dx} \left( AE \frac{du}{dx} \right) + 2x = 0 $$

This yields $$ 2x = 0 $$ which does not satisfy.

Checking the natural boundary condition, this equation is used

$$ \sigma \left(1 \right) = \left(E \frac{du}{dx} \right)_{x=1} = 0.1 $$

Substituting in $$ \alpha_1 $$ we get

$$ 0.1 - \frac{\frac{10}{3}}{A} = 0.1 $$

Thus, this boundary condition will satisfy only for large values of $$ A $$

=Problem 10: Fish and Belytschko Problem 3.4=

Repeat Problem 3.3 with the trial solution $$ u\left ( x \right ) = \alpha_{0} + \alpha_{1}\left ( x-3 \right ) + \alpha_{2}\left ( x-3 \right )^{2}$$.

Find
Consider a trial (candidate) solution of the form $$ u\left ( x \right ) = \alpha_{0} + \alpha_{1}\left ( x-3 \right ) + \alpha_{2}\left ( x-3 \right )^{2}$$ and a weight function of the same form. Obtain a solution to the weak form in Problem 3.1. Check the equilibrium equation in the strong form in Problem 3.1; is it satisfied?

Check the natural boundary condition; is it satisfied?

Given
The weak from is given by $$ \int_{1}^{3} \frac{\mathrm{d} w }{\mathrm{d} x}AE\frac{\mathrm{d} u}{\mathrm{d} x}dx = -0.1(wA)_{x=1} + \int_{1}^{3} 2xwdx, \forall w $$ with $$  w \left(3 \right) = 0 $$

Solution
The trial solution $$ u\left ( x \right )$$ and the weight function $$ w\left ( x \right )$$ can be written as follows:

$$ u\left ( x \right ) = \alpha_{0} + \alpha_{1}\left ( x-3 \right ) + \alpha_{2}\left ( x-3 \right )^{2}$$

$$ w\left ( x \right ) = \beta_{0} + \beta_{1}\left ( x-3 \right ) + \beta_{2}\left ( x-3 \right )^{2}$$

When:

$$ w \left(3 \right) = 0; \Rightarrow \beta_{0} = 0 $$

and for:

$$ u \left(3 \right) \Rightarrow \alpha_{0} = 0.001 $$

$$ u\left ( x \right )$$ and $$ w\left ( x \right )$$ can therefore be re-written as follows:

$$ u\left ( x \right ) = 0.001 + \alpha_{1}\left ( x-3 \right )+\alpha_{2}\left ( x-3 \right )^{2}$$

$$ w\left( x \right) = \beta_{1} \left(x - 3 \right)+\beta_{2}\left ( x-3 \right )^{2}$$

Whith their respective derivatives:

$$\frac{du}{dx} = \alpha_{1}+2\alpha_{2}\left(x-3\right)$$

$$\frac{dw}{dx} = \beta_{1}+2\beta_{2}\left(x-3\right)$$

Now substituting into the weak form yields:

$$ \int_{1}^{3} AE[\beta_{1}+2\beta_{2}\left(x-3\right)][\alpha_{1}+2\alpha_{2}\left(x-3\right)]dx = -0.1[A(\beta_{1}\left(x-3\right)+\beta_{2}\left(x-3\right)^{2})]_{x=1} + \int_{1}^{3} 2x[\beta_{1}\left(x-3\right)+\beta_{2}\left(x-3\right)^{2}]dx$$

which can be re-arranged as:

$$ \int_{1}^{3} AE[\beta_{1}+2\beta_{2}\left(x-3\right)][\alpha_{1}+2\alpha_{2}\left(x-3\right)]dx - 0.1[A(\beta_{1}\left(x-3\right)+\beta_{2}\left(x-3\right)^{2})]_{x=1} - \int_{1}^{3} 2x[\beta_{1}\left(x-3\right)+\beta_{2}\left(x-3\right)^{2}]dx = 0$$

Evaluating the integral we get:

$$ AE \left(2\alpha_{1}\beta_{1}-4\alpha_{1}\beta_{2}-4\alpha_{2}\beta_{1}+\frac{32\alpha_{2}\beta_{2}}{3}\right)- 0.1[A(4\beta_{2}-2\beta_{1})] - \left(8\beta_{2}-\frac {20\beta_{1}}{3}\right)=0 $$

Factoring out $$\left( \beta_{1}\right)$$ and $$\left( \beta_{2}\right)$$ we get:

$$\beta_{1}\left[2AE\left(\alpha_{1}-2\alpha_{2}\right)-0.2A+\frac{20}{3}\right]+\beta_{2}\left[4AE\left(\frac{8}{3}\alpha_{2}-\alpha_{1}\right)+0.4A-8\right]=0$$

$$\displaystyle {{\beta }_{1}}$$ and $$\displaystyle {{\beta }_{2}}$$ must be zero since they were arbitrarely chosen and because the above equation must hold. Therefore we can write both terms of the above equation equal to zeor as folows:

$$ \left[2AE\left(\alpha_{1}-2\alpha_{2}\right)-0.2A+\frac{20}{3}\right]=0$$

and

$$ \left[4AE\left(\frac{8}{3}\alpha_{2}-\alpha_{1}\right)+0.4A-8\right]=0 $$

Solving for $$\displaystyle {{\alpha }_{1}}$$ and $$\displaystyle {{\alpha }_{2}}$$ we find that:

$$\displaystyle {{\alpha }_{1}}=\frac{\frac{116}{3}+0.2A}{2AE} $$

and

$$\displaystyle {{\alpha }_{2}}=\frac{8}{AE} $$

The solution to the weak form can be written as


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

$$ u(x)= 0.001+\frac{\frac{116}{3}+0.2A}{2AE}(x-3)-\frac{8}{AE}(x-3)^{2} $$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * <p style="text-align:right;">$$\displaystyle (Eq. 10.1) $$
 * style= |
 * }

Checking if equilibrium in the strong form is satisfied, substituted into the equations:

$$ \frac{d}{dx} \left( AE \frac{du}{dx} \right) + 2x = 0 $$

Yields the following:

$$\displaystyle {2x-16=0}$$

which does not satisfy.

Checking if the natural boundary condition is satisfied, the following equation can be used:

$$ \sigma \left(1 \right) = \left(E \frac{du}{dx} \right)_{x=1} = 0.1 $$

=Problem 11: Solving using WRF=

Solving a differential equation using Weight Residual Form. From lecture 17

Given
The partial differential equation is:


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

$$ \frac{\partial}{\partial x}\left(2\frac{\partial u}{\partial x}\right) + 3 = 0 $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.1) $$
 * }
 * }

The boundary conditions are:

$$ u\left(1\right) = 0$$

$$\frac{du}{dx}\left(0\right) = -4 $$

As weighting funtion should be used:


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



$$ b_{i} =1+a_{i}\cos \left ( ix \right )+b_{i}\sin \left ( ix \right ) $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.2) $$
 * }
 * }

Find

 * 11.1 Find an approximate solution of the form $$U^{h}=\sum_{i=0}^{n}d_{i}b_{i}$$ for n=2
 * 11.2 Find two equations that enforce the boundary conditions
 * 11.3 Project the weight residues
 * 11.4 Display the equations in matrix form
 * 11.5 Solve for $$\displaystyle d $$
 * 11.6 Construct $$U^{h}$$ and plot $$U^{h}$$ and $$u$$
 * 11.7 Repeat 1.1 to 1.6 for n = 4 and n = 6

11.1 Find an approximate solution with n=2
When n = 2 we have:
 * {| style="width:100%" border="0"

$$ d_{i}=\begin{bmatrix} d0 & d1 & d2 & d3 & d4

\end{bmatrix} $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.3) $$
 * }
 * }

and


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

$$ b_{i}=\begin{bmatrix} 1 & \cos\left ( x\right) & \sin\left ( x\right) & \cos\left ( 2x\right) & \sin\left ( 2x\right)\end{bmatrix} $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.4) $$
 * }
 * }

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

$$ U^{h}=d0+d1\cos\left ( x \right )+d2\sin\left ( x\right )+ d3\cos\left ( 2x \right )+d4\sin\left ( 2x\right ) $$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * <p style="text-align:right;">$$\displaystyle (Eq. 11.5) $$
 * style= |
 * }

11.2 Find two equations that enforce the boundary conditions
The boundary conditions are:
 * {| style="width:100%" border="0"

$$ u\left(1\right)=0=d0+d1\cos\left(1\right)+d2\sin\left(1\right)+d3\cos\left(2\right)+d4\sin\left(2\right)$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * <p style="text-align:right;">$$\displaystyle (Eq. 11.6) $$
 * style= |
 * }

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

$$ \frac{du}{dx}(0)=-4 = d2+2d4$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * <p style="text-align:right;">$$\displaystyle (Eq. 11.7) $$
 * style= |
 * }

11.3 Project the weight residues
In this case changing the function u(x) by the approximated function $$U^{h}$$ in the PDE (Eq.11.1), we found that P(u) is:


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

$$ P\left ( U^{h} \right )=2\frac{\partial^2 U^{h}}{\partial x^2}+3=2d2+3 $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.8) $$
 * }
 * }

The (Eq.11.8) can be written as:


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

$$ P\left ( U^{h} \right )= -2\begin{bmatrix}0 & \cos\left(x\right) & \sin\left(x\right) & 4\cos\left(2x\right) &4\sin\left(2x\right)\end{bmatrix}\begin{bmatrix} d0\\ d1\\ d2\\ d3\\ d4\end{bmatrix}+3 $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.9) $$
 * }
 * }

Projecting the residue we have:


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

$$ \int_{a}^{b}b(x)P\left ( U^{h} \right )= 0 $$ After substitution of (Eq.11.4) and (Eq.11.9) in (Eq.11.10) the following equation is obtained
 * <p style="text-align:right;">$$\displaystyle (Eq. 11.10) $$
 * }
 * }
 * {| style="width:100%" border="0"

$$ 2\int_{0}^{1}\begin{bmatrix} 1\\ \cos\left (x\right ) \\ \sin\left ( x\right )\\ \cos\left (2x\right ) \\ \sin\left ( 2x\right ) \end{bmatrix}\begin{bmatrix} 0 &\cos\left(x\right) & \sin\left(x\right)&\cos\left(2x\right) & \sin\left(2x\right)\end{bmatrix}dx\begin{bmatrix} d0\\ d1\\ d2\\ d3\\ d4\end{bmatrix}=3\int_{0}^{1} \begin{bmatrix} 1\\ \cos\left ( x\right ) \\
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |

\sin\left ( x \right ) \\ \cos\left ( 2x\right ) \\

\sin\left ( 2x \right ) \end{bmatrix}dx $$


 * <p style="text-align:right;">$$\displaystyle (Eq. 11.11) $$
 * style= |
 * }

11.4 Display the equations in matrix form
Performing the product and then the integration indicated in (Eq.1.11) we have:
 * {| style="width:100%" border="0"

$$ \begin{bmatrix} 0 &1.68 & 0.91 &3.63& 5.66\\ 0 & 1.45& 0.70 & 3.55&4.49\\ 0 &0.70 & 0.54 &0.81 & 3.17\\ 0 & 0.88& 0.20 & 3.24& 1.65\\ 0 & 1.12& 0.79& 1.65& 4.75\end{bmatrix}\begin{bmatrix} d0\\ d1\\ d2\\ d3\\ d4\end{bmatrix}=\begin{bmatrix} 3\\ 2.52\\ 1.37\\ 1.36\\ 2.12\end{bmatrix} $$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * <p style="text-align:right;">$$\displaystyle (Eq. 11.12) $$
 * style= |
 * }

The latter equation is clearly of the form $$\displaystyle Kd=F$$. We can observe that K is not symmetric.

11.5 Solve for $$ d $$
In (Eq.11.12) we have three equations; but, as we need to enforce the boundary conditions we take the only one of them and solve it together with (Eq.11.6) and (Eq.11.7). From these equations we now have the system:


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

$$ \begin{bmatrix} 0 &\cos\left(1\right)   & \sin\left(1\right) &\cos\left(2\right)   & \sin\left(2\right)\\ 0 &0 &1 &0  & 2\\ 0  & -1.68  & -0.91 & -3.63  & -5.66\\ 0 &-1.45 & -0.70 & -3.55   & -4.49\\ 0 & -0.70& -0.54 & -0.81   & -3.17 \end{bmatrix}\begin{bmatrix} d0\\ d1\\ d2\\ d3\\

d4\end{bmatrix}=\begin{bmatrix} 0\\ -4\\ -3\\ -2.52\\ -1.37\end{bmatrix}$$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.13) $$
 * }
 * }

After solving the system in (Eq.11.13), by using the code shown in apedix, we obtain d:
 * {| style="width:100%" border="0"

$$\begin{bmatrix} d0\\ d1\\ d2\\ d3\\ d4\end{bmatrix}=\begin{bmatrix} 0.67\\ 4.88\\ -4.72\\ -0.80\\ -0.36\end{bmatrix}$$
 * style="width:2%; padding:1px; border:2px solid #FF0000" |
 * style="width:2%; padding:1px; border:2px solid #FF0000" |


 * <p style="text-align:right;">$$$$
 * style= |
 * }

11.6 Construct $$U^{h}$$ and plot $$U^{h}$$ and $$u$$
Replacing d in (Eq.11.5) we obtain the approximated solution:


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

$$ U^{h}=0.67+4\cos\left(x\right)-4.72\sin\left(x\right)-0.8\cos\left(2x\right)-0.36\sin\left(2x\right) $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.14) $$
 * }
 * }

On the other hand, the exact solution of the PDE (Eq.11.1) with the given boundary condition is:


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

$$ u\left ( x \right )=-\frac{3}{4}x^{2}-4x+\frac{19}{4} $$
 * <p style="text-align:right;">$$\displaystyle (Eq.11.15) $$
 * }
 * }

11.7 Repeat 11.1 to 11.6) for i = 4 and i = 6
Now we repeat the procedure using i = 4 and i=6 in Eq.11.4

As the procedure is the same shown before, we develop a routine using Matlab to perform the products, the integration and to solve the system of equations. The routine is show in the appendix.

In Fig.11.1 we show the exact and the approximated solution for i = 2,4, and 6. In that figure, it can be observed that all the numerical solutions fits very well the analytical solution.

We found a small error for i = 2 (0.0003); however the error lightly increases for i = 4 and decreases again for i = 6. As all the solutions are close to the analytical, we believe that this error is mainly due to round-off error which is increased with i = 4 and i = 6 due to the higher number of operations that the program has to perform. Fig. 11.2 shows the error variation as a function of i





Appendix
 Matlab Code: 

=Contributing Members=
 * James Roark
 * Abhijeet Bhalerao
 * Fernando Casanova
 * Yu Hao
 * Harrison Sheffield
 * Cedric Adam