text:
User:EML4500.f08.JAMAMA/FE/smith

ɥ$$ \ne \quad $$

http://meta.wikimedia.org/wiki/Help:Formula

=7.3: Newton's Law of Colling=

Given
Data Set G2DM2.0/D1

$$ \Omega = Domain \quad $$

$$ K=I, f(x,t)=0$$

$$ \frac {\partial u}{\partial t} = 0, ~g=2, ~h=3$$

$$ H=.5 ~u_\infty =.1$$

$$ \partial \Omega = \Gamma_g \bigcup \Gamma_h \bigcup \Gamma_H $$

and, The Partial Differential Equation is:

Find
1) The weak form for heat conduction in 2-D with boundary convection. 2) Develop semi-discrete equations (ODEs). Give detailed expression of all quantities. 3) Solve G2DM2.0/D1 using 2D LIBF till $$\displaystyle 10^{-6} $$ accuracy at the center. Plot solution in 3D with contour.

1.The Weak Form for heat conduction in 2D with boundary convection
From the PDE (Eqn 6.10.1) the weak form can be written as,

Using integration by parts,

Substituting into the integral we obtain,

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

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

\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{\frac{\partial }{\partial {{x}_{i}}}\left( w{{K}_{ij}}\frac{\partial u}{\partial {{x}_{j}}} \right)}d\text{ }\!\!\Omega\!\!\text{ }-\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{\frac{\partial w}{\partial {{x}_{i}}}\left( {{K}_{ij}}\frac{\partial u}{\partial {{x}_{j}}} \right)}d\text{ }\!\!\Omega\!\!\text{ +}\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{wf\left( x,t \right)}d\text{ }\!\!\Omega\!\!\text{ }-\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{w\rho c\frac{\partial u}{\partial t}}d\text{ }\!\!\Omega\!\!\text{ }=0 $$ $$
 *  $$ \displaystyle (Eq. 6.10.3)
 * }

Applying the Gauss Theorem on the first term gives,

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

\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}w{{K}_{ij}}\frac{\partial u}{\partial{{x}_{j}}}{{n}_{i}}d\text{ }\!\!\Gamma\!\!\text{ }-\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{\frac{\partial w}{\partial {{x}_{i}}}\left( {{K}_{ij}}\frac{\partial u}{\partial {{x}_{j}}} \right)}d\text{ }\!\!\Omega\!\!\text{ +}\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{wf\left( x,t \right)}d\text{ }\!\!\Omega\!\!\text{ }-\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{w\rho c\frac{\partial u}{\partial t}}d\text{ }\!\!\Omega\!\!\text{ }=0 $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle (Eq. 6.10.4)
 * }

The Heat Flux is given by,

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

q_{i}=-K_{ij}\frac{\partial u}{\partial x_{j}} $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle
 * }

Substituting the heat flux and rearranging the terms in Eq. 6.10.4 gives,

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

\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{-w{{q}_{i}}{n}_{i}}d\text{ }\!\!\Gamma\!\!\text{ }-\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{\frac{\partial w}{\partial {{x}_{i}}}\left( {{K}_{ij}}\frac{\partial u}{\partial {{x}_{j}}} \right)}d\text{ }\!\!\Omega\!\!\text{ +}\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{wf\left( x,t \right)}d\text{ }\!\!\Omega\!\!\text{ }-\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{w\rho c\frac{\partial u}{\partial t}}d\text{ }\!\!\Omega\!\!\text{ }=0 $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle (Eq. 6.10.5)
 * }

The first term can be rewritten as,

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

\int\limits_{\Gamma }{-w{{q}_{i}}{{n}_{i}}d\Gamma =}\int\limits_{-w{{q}_{i}}{{n}_{i}}d{{\Gamma }_{g}}}-\int\limits_{-w{{q}_{i}}{{n}_{i}}d{{\Gamma }_{h}}}-\int\limits_{-w{{q}_{i}}{{n}_{i}}d{{\Gamma }_{H}}} $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle
 * }

Select $$ \quad w $$ such that $$\quad w=0 $$ on $$\quad {{\Gamma }_{g}}$$

Organizing, Eq. 6.10.5 then becomes the following weak form

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

\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{w\rho c\frac{\partial u}{\partial t}d\text{ }\!\!\Omega\!\!\text{ +}}\int\limits_{\text{ }\!\!\Omega\!\!\text{ }}{\frac{\partial w}{\partial {{x}_{i}}}\left( {{K}_{ij}}\frac{\partial u}{\partial {{x}_{j}}} \right)d\text{ }\!\!\Omega\!\!\text{ }}+\int\limits_{wHu}d{{\text{ }\!\!\Gamma\!\!\text{ }}_{H}}=-\int\limits_{whd{{\text{ }\!\!\Gamma\!\!\text{ }}_{h}}}+\int\limits_{wH{{u}_{\infty }}d{{\text{ }\!\!\Gamma\!\!\text{ }}_{H}}}+\int\limits_{\Omega }{wf\left( x,t \right)d\text{ }\!\!\Omega\!\!\text{ }},$$ $$ $$ \quad \forall w $$ such that $$\displaystyle w=0 \quad on \quad {{\Gamma }_{g}} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle (Eq.6.10.6)
 * }

2) Develop Semidiscrete Equations (ODE's)
Let $$ u=[ N]{{u}_{e}}\quad $$ and $$ \quad w=[ N ]{{w}_{e}},$$ where $$\displaystyle [N] $$ is the Lagrange polynomial and $${u}_{e}~ and ~ {w}_{e} $$ are the nodal values of $$ u \quad $$ and $$w \quad $$ respectively. Now, take the derivative of basis functions with respect to x1 and x2

Substituting the Langrangian form into the weak form (Eq 6.10.6), yields the semi-dicrete form:

From Eq. 6.10.8, we can obtain the Heat Source vector, Conductivity matrix, and Capacitance Matrix.

The heat source vector is $$\tilde{f}(w)= -\int\limits_{h{{[N]}^{T}}d{{\Gamma }_{h}}}+\int\limits_{H{{[N]}^{T}}{{u}_{\infty }}d{{\Gamma }_{h}}}+\int\limits_{\Omega }{{{[N]}^{T}}f(x,t)d\Omega } $$

The Conducitivity matrix is $$\tilde{k}(w,u)= \int\limits_{\Omega }{{{[B]}^{T}}[B]d\Omega }+\int\limits_{H{{[N]}^{T}}[N]d{{\Gamma }_{h}}} $$

The Capacitance Matrix is $$\tilde{m}(w,u)= \int\limits_{\Omega }{\rho c{{[N]}^{T}}[N]\frac{\partial u}{\partial t}d\Omega } $$

3.Solve G2DM2.0/D1 Using 2D LIBF till $$10^{-6}$$ accuracy at center. Plot Solution in 3D w/Contour.
To solve the problem using LIBF, divide the x-axis and Y-axis into 2 parts 2D LIBF

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

$$ N_I (x,y) = L_{i,m}(x).L_{j,n}(y) \quad Where \ I = i+(j-1)m $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle                                                                     $$
 * }

Where,

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

$$ L_{i,m}(x) = \prod_{k=1\neq i}^{m}\frac{x-x_k}{x_i - x_k} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }
 * {| style="width:100%" border="0"

$$ L_{j,n}(x) = \prod_{k=1\neq j}^{m}\frac{y-y_k}{y_j - y_k} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

Consider a case of 9 nodes as shown below

The trial solution is

The trial solution must satisfy the essential boundary condition. Throughout boundary $$ \Gamma_g, \quad $$ the temperature should be a constant 2. Thus $$ d_1, d_2, d_3, d_6, d_9 =2. \quad $$ The matlab code used is shown and explained below

title: User:EML4500.f08.JAMAMA/FE/smith


text:
User:EML4500.f08.JAMAMA/FE/mcdaniel

= 7.3 =

Given/Find: Static Solution of the Circle

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

Let Omega a circle with a one-unit radius. Find the static solution such that: $$T=4$$ $$f(\underline x)=1$$ in Omega $$g=0$$ on $$\Gamma _g=\partial \Omega $$
 * style="width:95%" |
 * style="width:95%" |

all to 10^-6 accuracy at the center. Quadrilateral and Triangular elements should be used but for only 1/4 of the circle. Plot the deformed shape in 3D. $$
 *  $$ \displaystyle
 * }

Solution: Static Solution of the Circle Using Quadrilateral and Triangular Elements to an Accuracy of 10^-6

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

The Partial Differential Equation for a vibrating membrane is the following via transient analysis: $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle
 * }

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

T(\textrm{div \ grad \ u})+\textrm{f}=\rho \frac{\partial ^2u}{\partial t^2} $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle
 * }

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

Where: $$ u(\underline x,t)=$$ Transverse Displacement $$T=\quad$$ Tension (Force/Length) = Constant $$\textrm{f}(\underline x,t)=$$ Distributed Transverse Force $$\rho (\underline x)=$$ Mass Density (Mass/Area) & $$\textrm{div \ grad\ u}=\frac{\partial }{\partial x_i}\frac{\partial u}{\partial x_j}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle
 * }

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

The Weighted Residual Form is as follows: $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle
 * }

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

\int_{\Omega } w\left [ T\frac{\partial ^2u}{\partial x_i\partial x_j}+\textrm{f}-\rho \frac{\partial ^2u}{\partial t^2} \right ]d\Omega =0 $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle
 * }

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

Which can thereby be broken into three separate terms: $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle
 * }

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

\int_{\Omega } wT\frac{\partial ^2u}{\partial x_i\partial x_j}d\Omega =\mathbb{K} $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *  $$ \displaystyle
 * }

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

\int_{\Omega } w\textrm{f}\ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

\int_{\Omega } w\rho \frac{\partial ^2u}{\partial t^2}\ d\Omega=\tilde{m}(w,u) $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

This final equation will be used later in the Weak Form. For now, the first term above (set equal to $$\mathbb{K}$$), is equivalent to: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

\mathbb{K}=\int_{d\Omega } \frac{\partial }{\partial x_i}(wT\frac{\partial u}{\partial x_j})\ d\Omega-\int_{\Omega } \frac{\partial w}{\partial x_i}T\frac{\partial u}{\partial x_j}\ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

And by using divergence theorem, the first term of the equation above may be converted to yield: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

\mathbb{K}=\int_{\partial \Omega } n_iwT\frac{\partial u}{\partial x_j}d(\partial \Omega)-\int_{\Omega } \frac{\partial w}{\partial x_i}T\frac{\partial u}{\partial x_j}\ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

Now since $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

\partial \Omega =\Gamma _g\cup \Gamma _h\ \textrm{but }\ w=0\ \ \textrm{on}\ \ \Gamma _g $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

and the following are true: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

\frac{\partial w}{\partial x_i}=\textrm{grad}\ w=\bigtriangledown w $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

\frac{\partial u}{\partial x_j}=\textrm{grad}\ u=\bigtriangledown u $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

Then the K equation becomes: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

\mathbb{K}=\int_{\Gamma _h}(n_i\ wT\frac{\partial u}{\partial x_j}\ \partial \Gamma _h)-\int_{\Omega }\bigtriangledown wT\bigtriangledown u\ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

Otherwise written as: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

\mathbb{K}=\int_{\Gamma _h}w\ \underline nT\cdot \textrm{grad \ u}\ \partial \Gamma _h-T\int_{\Omega }\bigtriangledown w\bigtriangledown u\ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

where the last term is conveniently used to write ultimately the weak form of the weighted residual form which is the following: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

\tilde{\textrm{f}}(w)=\tilde{m}(w,u)+\tilde{k}(w,u)\ \ \ \forall w \ \ \ s.t.\ \ \ w=0\ \ \ on\ \ \ \Gamma _g $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

\tilde{m}(w,u)=\int_{\Omega } w\rho \frac{\partial ^2u}{\partial t^2}\ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

\tilde{k}(w,u)=T\int_{\Omega }\bigtriangledown w\bigtriangledown u \ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

\tilde{\textrm{f}}(w)=\int_{\Gamma _h}w \underline n \ T\ \textrm{grad u}\ \partial \Gamma _h+\int_{\Omega }w\textrm{f}\ d\Omega $$ $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 *
 $$ \displaystyle
 * }

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

From the Discrete Weighted Form, $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$\underline M_{FF}\underline d_F+\underline K_{FF}\underline d_F=\underline F_F-(\underline M_{FE}g+\underline K_{FE}g)$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

But because g=0 from the boundary conditions, as well as M due to the static requirement, the displacements may be derived: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$\underline K \underline d=\underline F$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$\therefore \underline d=\underline K^{-1}\underline F$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

Now to solve for the displacements: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$N_I(x,y)=L_{i,m}(x)\cdot L_{j,n}(y) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

Where $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$I=i+(j-1)m $$ and $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$L_{i,m}(x)=\prod_{K=1, K\neq i}^{m}\left ( \frac{x-x_k}{x_i-x_k} \right ) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

And setting the Shape Function to the following: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$N=\begin{bmatrix} \frac{1}{4}(1-\xi _1)(1-\xi _2)\\ \frac{1}{4}(1+\xi _1)(1-\xi _2)\\ \frac{1}{4}(1-\xi _1)(1+\xi _2)\\ \frac{1}{4}(1-\xi _1)(1+\xi _2) \end{bmatrix} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

The deformed shape coordinates may further be derived from these following two equations, where the summations are taken across all coordinate points taken from Abaqus output files: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$X_{1}^{e}=\varphi_{1}^{e}(\xi )=\sum_{I=1}^{nel}N_I(\xi )x_{1,I}^{e} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$(y=)X_{2}^{e}=\varphi_{2}^{e}(\xi )=\sum_{I=1}^{nel}N_I(\xi )x_{2,I}^{e} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

Using the following Jacobian matrix: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$J=\left [ \frac{\partial x_i}{\partial \xi _j} \right ]=\begin{bmatrix} \frac{\partial x_1}{\partial \xi _1} & \frac{\partial x_1}{\partial \xi _2}\\ \frac{\partial x_2}{\partial \xi _1} & \frac{\partial x_2}{\partial \xi _2}\\ \end{bmatrix} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

The following equations may now be utilized: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$\underline B_I(\xi )=\bigtriangledown _xN_{I}^{e}(\xi )=(\underline J^e)^{-T}(\xi )\bigtriangledown_\xi N_{I}^{e}(\xi ) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$k_{I3}^{e}=\int_{\omega =el}N_I(\xi )\cdot T\cdot \underline B_J (\xi )J(\xi )d(el) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$F_i=\int_{\omega =el}N_I(\xi )f\ det(J(\xi ))d(el) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

Implementing MATLAB at this point, in conjunction with output data files from Abaqus which contain element coordinate points for three separate quadrilateral meshes and three separate triangular meshes (included), the final displacements may be derived, to an accuracy of 10^-6: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

The above MATLAB script for the 6 different meshes generate the following results: $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

title: User:EML4500.f08.JAMAMA/FE/mcdaniel


text:
User:Eml5526.s11.team2.brenner/hwk2

=Problem 2.4: =

Solution
=Problem 2.6:=

Solution
=Problem 2.8: =

Given
Equations provided below can be found in Eqn (1) & (2) mtg. 10-4

Find
Show that 2.8.1 and 2.8.2 are equivalent.

Ensure that $$\displaystyle {u^h}(x) $$ satisfies the essential and natural boundary conditions.

where 2.8.4 is the essential boundary condition and 2.8.5 is the natural boundary condition

(Method I)
Taking the derivative with respect to x of both sides of (2.8.3) provides the following

Following dividing (2.8.6) through by $$\displaystyle P({u^h}(x)) $$ leads to

The weighted function $$\displaystyle {w^h}(x) $$ per Eqn (1) mtg. 10-3 takes the following form

Substituing (2.8.8) into (2.8.6) provides the following

Selecting a linear polynomial basis $$ \displaystyle \{ 1,x,{x^2}\} $$ on which to expand (2.8.9) and making use of the fact that (2.8.1) & (2.8.2) are equal to zero

Recall the requirement for the weighting function to be zero for $$\displaystyle w^h(x=1) = 0 $$ implies that

Making use of the Bubnov-Galerkin Method and we can construct a trial solution $$\displaystyle u^h(x) $$ of the form

Taking into account (2.8.8) and the requirement that the weighting function satisfy $$\displaystyle w^h(x=1) = 0 $$ reduces (2.8.12) to the following

Thereby satisfying (2.8.4) and the requirement on the essential boundary condition.

(Method B)
Begin by demonstrating that (8.2.1) is equivalent to (8.2.2). Assume the following finite approximation for the weighting function $$\displaystyle {w^h}(x) $$ which can be found in   Eqn (1) mtg. 10-3

Substituting (2.8.6) into (2.8.1) leads to the following expression

Expanding the right hand side of the summation found in (2.8.7) leads to the following form

(2.8.8) can be represented as a sum of integrals on $$\displaystyle \Omega $$

Next unexpand the summation of integrals in (2.8.9) and group the coefficients of $$\displaystyle \sum\limits_{i = 1}^n $$

Enforcing the elimination of the trivial solution requires that $$\displaystyle c_i \ne 0 $$. This leads to the following which satisfies (2.8.2).

Recall the requirement for the weighting function to be zero, which is as follows

Making use of the Bubnov-Galerkin Method we can construct a trial solution $$\displaystyle u^h(x) $$ of the form

Taking into account (2.8.13) and the requirement that the weighting function satisfy $$\displaystyle w^h(x=1) = 0 $$ (2.8.13) reduces

Thereby satisfying (2.8.4) which is the requirement of the essential boundary condition.

Continue by demonstrating that (8.2.2) is equivalent to (8.2.1).

Multiple (8.2.2) by $$\displaystyle \sum\limits_{i = 1}^n $$ which provides the following

Reorganize (8.2.13) as follows

Per (2.8.6) substitute (2.8.16 ) into (2.8.15) which results in the following

Expand (2.8.17) to a difference of integrals $$\displaystyle \int_{\Omega} $$

Using (2.8.8) and (2.8.21) $$\displaystyle g(x) $$ can also be represented by

Substituting (2.8.19) into (2.8.18) provides the following

Assume the following finite approximation for the trial solution $$\displaystyle {w^h}(x) $$ which can be found in   Eqn (2) mtg. 10-3

Which modifies (2.8.19) into the following

substituting (2.8.6) into (2.8.22) satifies the proof.

title: User:Eml5526.s11.team2.brenner/hwk2


text: Category:Reviews

title: Category:Reviews


text: File:Nm1.s11.team5.HW2.P2.fig1.jpg .png

title: File:Nm1.s11.team5.HW2.P2.fig1.jpg .png


text:
User:Egm6341.spring-11.5.D/HW2

= Problem 1 - Taylor Series = = Problem 2 -Taylor Series =

Problem statement
Construct the first 10 Taylor series expansions of: $$ f(x)=sin x, \ x \in [0,\pi] $$

Solution
 Solved my self using Mathematica 7 and the following commands: Table[Series[Sin[x], {x, 3*Pi/8, n}], {n, 0, 10}]; MatrixForm[Normal[%]]

$$ \left( \begin{array}{c} \cos \left(\frac{\pi }{8}\right) \\ \left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ \left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ -\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ -\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)+\frac{1}{24} \left(x-\frac{3 \pi }{8}\right)^4 \cos \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ \frac{1}{120} \left(x-\frac{3 \pi }{8}\right)^5 \sin \left(\frac{\pi }{8}\right)-\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)+\frac{1}{24} \left(x-\frac{3 \pi }{8}\right)^4 \cos \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ \frac{1}{120} \left(x-\frac{3 \pi }{8}\right)^5 \sin \left(\frac{\pi }{8}\right)-\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)-\frac{1}{720} \left(x-\frac{3 \pi }{8}\right)^6 \cos \left(\frac{\pi }{8}\right)+\frac{1}{24} \left(x-\frac{3 \pi }{8}\right)^4 \cos \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ -\frac{\left(x-\frac{3 \pi }{8}\right)^7 \sin \left(\frac{\pi }{8}\right)}{5040}+\frac{1}{120} \left(x-\frac{3 \pi }{8}\right)^5 \sin \left(\frac{\pi }{8}\right)-\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)-\frac{1}{720} \left(x-\frac{3 \pi }{8}\right)^6 \cos \left(\frac{\pi }{8}\right)+\frac{1}{24} \left(x-\frac{3 \pi }{8}\right)^4 \cos \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ -\frac{\left(x-\frac{3 \pi }{8}\right)^7 \sin \left(\frac{\pi }{8}\right)}{5040}+\frac{1}{120} \left(x-\frac{3 \pi }{8}\right)^5 \sin \left(\frac{\pi }{8}\right)-\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)+\frac{\left(x-\frac{3 \pi }{8}\right)^8 \cos \left(\frac{\pi }{8}\right)}{40320}-\frac{1}{720} \left(x-\frac{3 \pi }{8}\right)^6 \cos \left(\frac{\pi }{8}\right)+\frac{1}{24} \left(x-\frac{3 \pi }{8}\right)^4 \cos \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ \frac{\left(x-\frac{3 \pi }{8}\right)^9 \sin \left(\frac{\pi }{8}\right)}{362880}-\frac{\left(x-\frac{3 \pi }{8}\right)^7 \sin \left(\frac{\pi }{8}\right)}{5040}+\frac{1}{120} \left(x-\frac{3 \pi }{8}\right)^5 \sin \left(\frac{\pi }{8}\right)-\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)+\frac{\left(x-\frac{3 \pi }{8}\right)^8 \cos \left(\frac{\pi }{8}\right)}{40320}-\frac{1}{720} \left(x-\frac{3 \pi }{8}\right)^6 \cos \left(\frac{\pi }{8}\right)+\frac{1}{24} \left(x-\frac{3 \pi }{8}\right)^4 \cos \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \\ \frac{\left(x-\frac{3 \pi }{8}\right)^9 \sin \left(\frac{\pi }{8}\right)}{362880}-\frac{\left(x-\frac{3 \pi }{8}\right)^7 \sin \left(\frac{\pi }{8}\right)}{5040}+\frac{1}{120} \left(x-\frac{3 \pi }{8}\right)^5 \sin \left(\frac{\pi }{8}\right)-\frac{1}{6} \left(x-\frac{3 \pi }{8}\right)^3 \sin \left(\frac{\pi }{8}\right)+\left(x-\frac{3 \pi }{8}\right) \sin \left(\frac{\pi }{8}\right)-\frac{\left(x-\frac{3 \pi }{8}\right)^{10} \cos \left(\frac{\pi }{8}\right)}{3628800}+\frac{\left(x-\frac{3 \pi }{8}\right)^8 \cos \left(\frac{\pi }{8}\right)}{40320}-\frac{1}{720} \left(x-\frac{3 \pi }{8}\right)^6 \cos \left(\frac{\pi }{8}\right)+\frac{1}{24} \left(x-\frac{3 \pi }{8}\right)^4 \cos \left(\frac{\pi }{8}\right)-\frac{1}{2} \left(x-\frac{3 \pi }{8}\right)^2 \cos \left(\frac{\pi }{8}\right)+\cos \left(\frac{\pi }{8}\right) \end{array} \right) /,$$

Residual
$$Max\left|R(x)_{n+1}\right|=\frac{\left(x-x_0\right){}^{n+1} Max\left|f(\zeta )^{n+1}\right|}{(n+1)!}, \zeta \in \left[x_o, x\right],x_0=\frac{3 \pi }{8},x=\frac{3 \pi }{4}$$ When $$n+1$$ is even $$ Max\left|f(\zeta )^{n+1}\right|=1$$ When $$n+1$$ is odd $$ Max\left|f(\zeta )^{n+1}\right|=\frac{1}{\sqrt{2}}$$ $$\therefore$$ When $$n+1$$ is even: $$Max\left|R(x)_{n+1}\right|=\frac{1.1781^{n+1}}{(n+1)!}$$ And When $$n+1$$ is odd: $$Max\left|R(x)_{n+1}\right|=\frac{1.1781^{n+1}}{\sqrt{2}(n+1)!}$$

Plot expansions
Homework 2 Problem 2 Figure 1: First 10 Taylor series expansions of Sin(x), from 0 to 2$$\pi$$

= Problem 3 - = = Problem 4 - = = Problem 5 - = = Problem 6 - = = Problem 7 - = = Problem 8 - =

= Problem 9 - =

= Problem 10 - =

= Problem 11 - = = Problem 12 - =

Function to gernate Lagrange Interpolation funcitons

LagrangeInterp[data_, x_] /; MatrixQ[data] := Flatten[Module[{n = Length[data], xl},{xl = data};Table[Product[(x - xli)/(xlj - xli), {i,Complement[Range[n], {j}]}], {j, n}]]\[Transpose]]

f[x_] := (Exp[x] - 1)/x

Mathematica code to generate evaluation points, data points, Interpolation funcitons,and interpolated functions for n=[1,2,4,8,16]

fn = Simplify[Table[fpts = Table[f[i], {i, -1.001, 1, 2.001/n}, {j, 1}]; xloc = Table[i, {i, -1.001, 1, 2.001/n}, {j, 1}]; LIf = LagrangeInterp[xloc, x]; LIf.fpts, {n, {1, 2, 4, 8, 16}}]]

Resulting Interpolation functions:

$$ \left( \begin{array}{c} 0.542941 x+1.17534 \\ 0.175144 x^2+0.543116 x+1.00002 \\ 0.00858167 \left(x^2-0.403782 x+11.4075\right) \left(x^2+5.46541 x+10.2149\right) \\ 2.8018942259677715*10^{-6} \left(x^2-5.53948 x+34.1119\right) \left(x^2+0.275742 x+21.7425\right) \left(x^2+5.34705 x+21.6285\right) \left(x^2+8.95483 x+22.2488\right) \\ 1.1816155165433884*10^{-8} x^{16}-1.6916601452976465*10^{-10} x^{15}-3.818422555923462*10^{-8} x^{14}-1.5643308870494366*10^{-10} x^{13}+4.784669727087021*10^{-8} x^{12}+1.760781742632389*10^{-9} x^{11}-6.548361852765083*10^{-9} x^{10}+2.74396370514296*10^{-7} x^9+2.7657661121338606*10^{-6} x^8+0.0000248015 x^7+0.000198411 x^6+0.00138889 x^5+0.00833333 x^4+0.0416667 x^3+0.166667 x^2+0.5 x+1. \end{array} \right)

$$

The difference between the integration of the original function and the interpolation funcitons is used to compare their reletive accurcy. The Mathematica code to accomplish this is:

N[Integrate[f[x], {x, -1, 1}]] - Integrate[fn, {x, -1, 1}]

The results are as follows:

$$\left( \begin{array}{c} -0.236179 \\ -0.00230395 \\ -9.742760151087992*10^{-6} \\ -1.1224932094933138*10^{-10} \\ 3.831686079536212*10^{-10} \end{array} \right)$$

= Problem 13 - =

title: User:Egm6341.spring-11.5.D/HW2


text:
User:EML4500.f08.JAMAMA/NM

= Cauchy Distribution =

Given
Cauchy Random Noise
 * {| style="width:100%" border="0"

$$\underline{F} = \underline{I} + \triangle \underline{A} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 7.0.1)
 * }

$$\underline{I} = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix},\quad \underline{A} = \begin{pmatrix} -0.2& 1 \\ -1 & -0.2 \end{pmatrix}, \quad \triangle = 0.02 \quad \underline{x_0}= \dbinom{{{x}_{0}}^{1}}{{{x}_{0}}^{2}} $$

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

$$LSSM: \quad \underline{x}_{k+1} = \underline{F} \quad \underline{x}_{k}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 7.0.2)
 * }

Find
Plot $$ ( \underline{{x}_{j}}, \quad j=0,1,2,...)$$ in single-slit diffraction

Solution

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

$$ LSSMRN: \quad \underline{{x}_{k+1}} = \underline{F} \quad \underline{x_k} + \underline{G} \quad \underline{w}_{k+1}$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 7.0.1)
 * }

Result is

Given
An ideal mass-spring-damper system with mass m, spring constant k and viscous damper of damping coefficient c is subject to a force F=u.
 * {| style="width:100%" border="0"

$$ \left\{ {\begin{array}{*{20}{c}} \dot x_1 \\ \dot x_2 \\ \end{array}} \right\} = \left[ {\begin{array}{*{20}{c}} 0 & 1 \\  -\frac{k}{m} & -\frac{c}{m} \\ \end{array}} \right] \left\{ {\begin{array}{*{20}{c}} x_1 \\ x_2 \\ \end{array}} \right\} + \left[ {\begin{array}{*{20}{c}} 0 & 0 \\  1 & 0 \\ \end{array}} \right] \left\{ {\begin{array}{*{20}{c}} u \\ 0 \\ \end{array}} \right\} $$   $$ (Eq.7.0.2) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

Finding
Let $$k=1,/ m=1/2,/ x_0=[0.8,-0.4]^T$$ For $$u=0.5$$ Cauchy noise, and $$c=\frac{3}{2}c_{cr}$$, plot $$x_k$$

Solution
Egm6341.s11.team3.sungsikkim 03:29, 21 April 2011 (UTC)

= =

Given: Verhulst Equation with initial conditions
Verhulst equation:

$$ \frac{d}{dt}x(t) = rx(1 - \frac{x}{x_{max}}) \ $$

$$ x_{max} = 15 \ $$ and $$ r = 1.4 \ $$

with initial condition of $$ x_0 \ $$

for time $$ t \in [0, 10] \ $$

Case 1:

$$ x_0 = 3 < \frac{1}{2} x_{max} \ $$

Case 2:

$$ x_0 = 9 > \frac{1}{2} x_{max} \ $$

Find: By using Hermite - Simpson rule, integrate given Verhulst equation
From given Prof. Vu Quoc's notes for

$$ z_{i+ 1} = z_i + \frac{h/2}{3} \left ( f_i + 4f_{i+\frac{1}{2}}[g(z_i,z_{i+1})],t_{i+\frac{1}{2}} + f_{i+1} \right ) \ $$

integrate this function

Solution
Main reference pages were meeting 38 and 39 notes

From given Prof. Vu Quoc's notes for

$$ z_{i+ 1} = z_i + \frac{h/2}{3} \left ( f_i + 4f_{i+\frac{1}{2}}[g(z_i,z_{i+1})],t_{i+\frac{1}{2}} + f_{i+1} \right ) \ $$

,where

$$ f_i = f(z_i, f_i) \ $$

$$ f_{i+1} = f(z_{i+1}, t_{i+1}) \ $$

$$ f_{i+\frac{1}{2}} = f[g(z_i,z_{i+1},t_{i+\frac{1}{2} })] = f(z_{i+\frac{1}{2}}, t_{i+\frac{1}{2}}) \ $$

,where $$ t_{i+\frac{1}{2}} = t_i + \frac{h}{2} \ $$ and as it was proven in HW 6.6 as well as given by Prof. Vu Quoc's

$$ z_{i+ \frac{1}{2}} = \frac{1}{2}(z_i + z_{i+1}) + \frac{h}{8}(f_i - f_{i+1}) \ $$

,note that $$ h \ $$ is the time step

Also specifying the equation to be of a form $$ F(z_{i+ 1})= 0 \ $$

$$ F(z_{i+ 1})= -z_{i+ 1} + z_i + \frac{h/2}{3} \left ( f_i + 4f([g(z_i,z_{i+1})],t_{i+\frac{1}{2}}) + f_{i+1} \right ) = 0\ $$

and, finally, as mentioned in

$$ z_{i+1}^{(k+1)} = z_{i+1}^{(k)} - \left [ \frac{d}{dz}F(z_{i+1}^{(k)}) \right ]^{-1} F(z_{i+1}^{(k)}) \ $$

in our case, for a given function, we have first order derivative. Therefore

$$ z_{i+1}^{(1)} = z_{i+1} - \left [ \frac{d}{dz}F(z_{i+1}) \right ]^{-1} F(z_{i+1}) \ $$

Then check the condition for convergence as it is shown in

$$ z_{i+1}^{(1)} - z_{i+1} ?< Tolerance \ $$

If condition is not met, choose $$ z_{i+1}^{(1)} = z_{i+2} \ $$ and iterate again to get $$ z_{i+2}^{(1)} \ $$ until the tolerance is met

Egm6341.s11.team3.rakesh 02:24, 21 April 2011 (UTC)

Egm6341.s11.team3.JA 17:00, 21 April 2011 (UTC)

Solution Part 2:
All of the necessary information is given in the problem statement. Therefore, all that is left to do is perform the calculations using the proposed methods. The same MatLAB code from Solution Part 1 may be used, with the exception of the file main.m. This file will be replaced by mainB.m as shown below, where the only difference between these files is that mainB.m uses the angle of attack governed by Equation (7.2.2.B), and mainB does not calculate J or ode45.

Below are the new resulting state variable plots, followed by the new M-File that produced them:

Egm6341.s11.team3.russo 03:50, 20 April 2011 (UTC)

title: User:EML4500.f08.JAMAMA/NM


text:
User:EML4500.f08.JAMAMA/NM/rakesh

rakesh

title: User:EML4500.f08.JAMAMA/NM/rakesh


text:
User:EML4500.f08.JAMAMA/NM/sungsikkim

sungsikkim

title: User:EML4500.f08.JAMAMA/NM/sungsikkim


text:
User:EML4500.f08.JAMAMA/NM/ren

=HW 5.4=

Given: An ideal mass-spring-damper system
An ideal mass-spring-damper system with mass m, spring constant k and viscous damper of damping coefficient c is subject to a force F=u.

Find: Obtain the governing equation and plot system respond
1) Derive equation of motion in terms of $$d,c,k,m,u$$.

2) Let $$ x=\left\{ {\begin{array}{*{20}{c}} d \\  \dot d  \\ \end{array}} \right\} =\left\{ {\begin{array}{*{20}{c}} x_1  \\  x_2  \\ \end{array}} \right\}  \ $$, find $$\left(F,G\right)$$

3) Find $$c_{cr}$$ in terms of $$k,m$$ such that system is critically damped.

4) Let $$k=1,/ m=1/2,/ x_0=[0.8,-0.4]^T$$

a) For $$u=0$$, plot $$x_k$$ for $$c=\frac{1}{2}c_{cr}, c_{cr}, \frac{3}{2}c_{cr}$$.  b) For $$u=0.5$$ Gaussian noise, and $$c=\frac{3}{2}c_{cr}$$, plot $$x_k$$ c) For $$u=0.5$$ Cauchy noise, and $$c=\frac{3}{2}c_{cr}$$, plot $$x_k$$

Solution
1)

For the spring k
 * {| style="width:100%" border="0"

$$ f_k=kd \ $$ $$ (Eq.5.4.1) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

For the damping c
 * {| style="width:100%" border="0"

$$ f_c=c\cdot \frac{d}{dt}d(t)=c\dot d \ $$ $$ (Eq.5.4.2) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

For the mass m
 * {| style="width:100%" border="0"

$$ m\frac{d^2}{dt^2}d(t) =u-f_k-f_c  \ $$ $$ (Eq.5.4.3) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

Combine (Eq.5.4.1-3), we obtain the governing equation of the system

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

$$ m\ddot d = u - c\dot d-kd \quad \Rightarrow\quad m\ddot d + c\dot d+kd= u  \ $$ $$ (Eq.5.4.4) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

2)

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

$$ \underline{x}=\left\{ {\begin{array}{*{20}{c}} d \\  \dot d  \\ \end{array}} \right\} =\left\{ {\begin{array}{*{20}{c}} x_1  \\  x_2  \\ \end{array}} \right\}  \ $$ $$ (Eq.5.4.5) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

From (Eq.5.4.4), we have

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

$$ \left\{ {\begin{array}{*{20}{c}} \dot x_1 = x_2 \\  \dot x_2 = u - \frac{c}{m}x_2-\frac{k}{m}x_1\\ \end{array}} \right. $$   $$ (Eq.5.4.6) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

that is

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

$$ \left\{ {\begin{array}{*{20}{c}} \dot x_1 \\ \dot x_2 \\ \end{array}} \right\} = \left[ {\begin{array}{*{20}{c}} 0 & 1 \\  -\frac{k}{m} & -\frac{c}{m} \\ \end{array}} \right] \left\{ {\begin{array}{*{20}{c}} x_1 \\ x_2 \\ \end{array}} \right\} + \left[ {\begin{array}{*{20}{c}} 0 & 0 \\  1 & 0 \\ \end{array}} \right] \left\{ {\begin{array}{*{20}{c}} u \\ 0 \\ \end{array}} \right\} $$   $$ (Eq.5.4.7) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

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

$$ \Rightarrow \quad \underline{A} = \left[ {\begin{array}{*{20}{c}} 0 & 1 \\  -\frac{k}{m} & -\frac{c}{m} \\ \end{array}} \right] \quad \underline{B}= \left[ {\begin{array}{*{20}{c}} 0 & 0 \\  1 & 0 \\ \end{array}} \right]$$ $$ (Eq.5.4.8) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 * }

According to the lecture notes, we have

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

$$ \underline{F}=\underline{I}+\Delta t \underline{A} = \left[ {\begin{array}{*{20}{c}} 1 & \Delta t \\ -\frac{k\Delta t}{m} & 1-\frac{c\Delta t}{m} \\ \end{array}} \right] \quad \underline{G}=\Delta t \underline{B}= \left[ {\begin{array}{*{20}{c}} 0 & 0 \\  \Delta t & 0 \\ \end{array}} \right]$$ $$ (Eq.5.4.9) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

3) For critical damping system, the governing equation can be written as
 * {| style="width:100%" border="0"

$$ \ddot d + 2\omega_0\dot d+\omega^2 d=const. \ $$   $$ (Eq.5.4.10) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

Comparing (Eq.5.4.10) with (Eq.5.4.4), we have

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

$$ \frac{k}{m} = \left(\frac{c_{cr}}{2m}\right)^2 \quad\Rightarrow\quad c_{cr}=2\sqrt{mk}\ $$ $$ (Eq.5.4.11) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

4) The governing equation (Eq.5.4.7) can be integrated according to

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

$$ \underline{x}_{k+1}=\underline{F}\underline{x}_k+\underline{G}\underline{u}_{k+1} \quad\ $$ $$ (Eq.5.4.12) \ $$ with  $$\quad x_0=[0.8,-0.4]^T \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

Matlab code is written to solve the problem and attached below:

Find: understand Kessler's code line by line
Using $$(p_{2i},\ p_{2i+1}),\ i=1,2,3$$ to understand Kessler's code line by line, starting with the best of Spring 2010

Solution
The Kessler's code is to compute the coefficients $$p_{2n}$$ and the constants $$c_{n}$$ where $$n=1,2,3,...,n$$ ,using the following algorithm:

1) $$i=0$$:
 * {| style="width:100%" border="0"

$$p_1(t)=c_1t,\ c_1=-1 $$ $$ (Eq.6.6.1) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

2) $$i=1,2,3,...$$:
 * {| style="width:100%" border="0"

$$c_{2i+1}=-\sum^{i-1}_{j=0}\frac{c_{2j+1}}{[2(i-j)+1]!}$$ $$ (Eq.6.6.2) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

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

$$p_{2i}(t)=\sum^{i}_{j=0}c_{2j+1}\frac{t^{2(i-j)}}{[2(i-j)]!}$$ $$ (Eq.6.6.3) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }
 * function traperror

The input of this function $$n$$ is a integer that ranges from 1 to $$\infty$$. The function returns two arrays $$c$$ and $$p$$, corresponding to the constants $$c_{n}$$ and the coefficients $$p_{2n}$$.

Set initial values for f, g, cn and cd, where the initial value of cn and cd are the numerator and denominator of $$c_1=-1$$, respectively.

Start iteration from $$k=1$$ to $$k=n$$

Update vector f according to $$f^k_i=f^{k-1}_ig^{k-1}_i(g^{k-1}_i + 1)$$, where $$f^k_i$$ is the denominator of the factors of $$c_i$$ in (Eqn 6.6.2) as $$i=k$$

Calculate $$c_{2k-1}$$ according to (Eqn 6.6.2)

Update vector f and g, where $$g_i$$ and $$f_i$$ are respectively the numerator and denominator of the factors of $$c_{2i+1}$$ in (Eqn 6.6.2) as $$i=k$$

Calculate $$p_{2k}$$ according to (Eqn 6.6.3)

Add $$c_{2k+1}$$ and $$p_{2k}$$ to the output arrays c and p

 * function fracsum

This function is to calculate the summation of $$n_i/d_i$$ from vectors $$n$$ and $$d$$. The outputs nsum and dsum returns the numerator and denominator,respectively. $$\frac{nsum}{dsum}=\sum_{i}\frac{n_i}{d_i}$$

The function round is used to round off a number to its nearest integer. The function gcd returns an array containing the greatest common divisors of the corresponding elements of the two input arrays.

This is to remove the common divisors from $$n_i$$ and $$d_i$$, while the ratio $$n_i/d_i$$ is preserved.

The function lcm returns the least common multiple of corresponding elements of two input arrays

Given: inconsistent Trapezoidal Simpson algorithm

 * inconsistent Trapezoidal Simpson algorithm

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

$$ z_{i+\frac{1}{2}}=\frac{1}{2}(z_{i}+z_{i+1}) \ $$ $$ (Eq.7.3.1) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

 * Hermite - Simpson algorithm

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

$$ z_{i+\frac{1}{2}}=\frac{1}{2}(z_{i}+z_{i+1})+\frac{h}{8}(f_{i}-f_{i+1}) \ $$ $$ (Eq.7.3.2) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

Find: Solve the logistic equation

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

$$ \frac{dz}{dt}=rx(1-\frac{x}{x_{max}}) \ $$ $$ (Eq.7.3.3) \ $$
 * style="width:95%" |
 * style="width:95%" |
 *

 * }

with $$x_{max}=15$$, $$r=1.4$$, and two initial conditions

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

$$x_0=3<\frac{1}{2}x_{max}$$
 * style="width:95%" |
 * style="width:95%" |

$$x_0=9<\frac{1}{2}x_{max}$$

$$t\in [0,20]$$
 * }
 * }

Solution
Then complete Hermite-Simpson time-stepping algorithm is given in. The only difference of the inconsistent Trapezoidal Simpson algorithm from the Hermite-Simpson algorithm is that the value of $$z_{i+1/2}$$ is evaluated as the mean value of $$z_{i}$$ and $$z_{i+1}$$, as shown in (Eq.7.3.2).

From observation of the above figures, the solutions given by the two algorithms almost coincide, regardless of the different initial conditions. Yet, the Hermite-Simpson algorithm is superior in the sense that it gives better estimation to the function value at the middle point.

128.227.113.35 19:25, 21 April 2011 (UTC)

title: User:EML4500.f08.JAMAMA/NM/ren


text:
User:EML4500.f08.JAMAMA/NM/russo

russo

title: User:EML4500.f08.JAMAMA/NM/russo


text:
User:EML4500.f08.JAMAMA/NM/Xia

= Numerical Methods =

Extra Solution
From Taylor series :

$$ e^x = 1 + x + \frac{x^2}{2} + \frac{x^3}{6} + \frac{x^4}{24} + \frac{x^5}{120} + ... $$

Then the original fuction could be expressed as:

$$ f(x) = 1 + \frac{x}{2} + \frac{x^2}{6} + \frac{x^3}{24} + \frac{x^4}{120} + ... $$

$$ \Rightarrow \lim_{x\to 0} f(x) = 1 + \frac{0}{2} + \frac{0^2}{6} + \frac{0^3}{24} + \frac{0^4}{120} + ... =1 $$

The original function is plotted from 0 to 1. For comparison, serie-form of f(x) of order 2nd and 3rd are also shown in the following plot.

Extra Solution
From Taylor series :

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

$$ e^x = 1 + \sum_{i=1}^n \frac{x^i}{i!} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle
 * }

Then the original fuction could be expressed as:

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

$$ f(x) = \frac{e^{x}-1}{x} = \sum_{i=1}^n \frac{x^{i-1}}{i!} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle
 * }

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

$$ \Rightarrow P_n(x) = \sum_{i=0}^{n} \frac{x^{i}}{(i+1)!} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 2.12)
 * }

To find $$ R_n(x) $$, we need to solve for $$ f^{(n+1)}(x) $$ which is difficult. So, here we compute $$ R_n(x) $$ from the previous function $$ e^{x} $$. From (Eq. 2.7), we could write $$ e^{x} $$ as:

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

$$ e^{x} = 1+ \sum_{i=1}^{n+1} \frac{x^i}{i!} + \frac{x^{n+2}}{(n+2)!} e^\xi, \quad \xi \in [0, x] $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle $$
 * }

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

$$ \Rightarrow f(x) = \frac{e^{x}-1}{x} = \sum_{i=0}^{n} \frac{x^i}{(i+1)!} + \frac{x^{n+1}}{(n+2)!} e^\xi, \quad \xi \in [0, x] $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 2.13)

$$
 * }

Comparing (Eq. 2.12) and (Eq. 2.13) yields,

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

$$ R_n(x) = \frac{x^{n+1}}{(n+2)!} e^\xi, \quad \xi \in [0, x] $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle

$$
 * }

Solution
Noting from equation 3.10

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

$$ f_n(x) = {P}_{n}(x) = \sum_{i=0}^{n} \frac{x^{i}}{(i+1)!} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 4.1)
 * }

So the integration could be derived as

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

$$ I_n = \int_{-1}^{1} f_n(x)dx = \int_{-1}^{1} \sum_{i=0}^{n} \frac{x^{i}}{(i+1)!}dx = \sum_{i=0}^{n} \int_{-1}^{1}\frac{x^{i}}{(i+1)!}dx $$ $$ substitute $$ \quad n $$ with $$ \quad 2k $$ yields,
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 4.2)
 * }

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

$$ I_n = \sum_{i=0}^{k} \frac{2}{(2i+1)(2i+1)!} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.3)
 * }

As $$ \quad n = 2,4,8,16,... $$, use matlab to solve for the integration and the code is attached below:

The results are shown in the following table:

Solution
From Simpson's rule :

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

$$ I_n(x) = \frac {h}{3} [f_0 + 4f_1 + 2\sum_{i=2}^{n-2} f_i + 4f_{n-1} + f_n] $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

Matlab code is written to solve the problem and attached below:

The solution and iteration history is shown in following table:

Solution
For Gauss-Legendre quad, the integration could be written as :

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

$$ I_n(x) = \sum_{i=1}^{n} w_if(x_i) $$ for $$ i=1,2,...,n $$ and $$ x_i $$ are the roots for (Eq) $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

The weight functions are given by :

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

$$ w_i = \frac {-2}{(n+1)P_{n}^{(1)}(x_i)P_{n+1}(x_i)} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ P_{n}(x) = \sum_{i=0}^{[n/2]} (-1)^i\frac{(2n-2i)!x^{n-2i}}{2^ni!(n-i)!(n-2i)!} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

Combining the above equations, we can obtain the roots and weight function values for n=2 and n=4 as presented in lecture notes. The following matlab code uses those values directly and compute the integrations for n=2 and n=4.

The results are shown in the following table:

Obviously from the above table, the error for n=4 already smaller than the order of 0.000001, so there is no need of further computaions.

Problem
Given from lecture notes:

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

$$ p_2(x)=\sum_{i=0}^{n=2} l_{i,2}(x)f(x_i) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.1)
 * }
 * {| style="width:100%" border="0"

$$ l_{i,n}=\prod_{j=0,j \neq i}^{n} \frac{x-x_j}{x_i-x_j} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.2)
 * }

Derive: Simple simpson's rule.

Solution
From definition:

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

$$ I_2(f) = \int_{a}^{b} f_2(x)dx = \sum_{i=0}^{2} \int_{a}^{b} l_{i,2}(x)dx f(x_i) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.3)
 * }

With

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

$$ x_0=a;x_1=\frac {a+b}{2};x_2=b. $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

from (Eq. 11.2), we can obtain:

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

$$ l_{0,2}(x) = \frac {(x-x_1)(x-x_2)}{(x_0-x_1)(x_0-x_2)} = \frac {2}{(a-b)^2} [x^2-\frac{a+3b}{2} x +\frac{(a+b)b}{2}] $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.4)
 * }

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

$$ l_{1,2}(x) = \frac {(x-x_0)(x-x_2)}{(x_1-x_0)(x_1-x_2)} = -\frac {4}{(a-b)^2} [x^2-(a+b) x +ab] $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.5)
 * }

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

$$ l_{2,2}(x) = \frac {(x-x_0)(x-x_1)}{(x_2-x_0)(x_2-x_1)} = \frac {2}{(b-a)^2} [x^2-\frac{3a+b}{2} x +\frac{(a+b)a}{2}] $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.6)
 * }

Integrate (Eq. 11.4)(Eq. 11.5)(Eq. 11.6), we can solve

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

$$ \int_{a}^{b}l_{0,2}(x) = \frac {2}{(b-a)^2} [-\frac{1}{12} (a-b)^3] = \frac{1}{6}(b-a) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.7)
 * }

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

$$ \int_{a}^{b}l_{1,2}(x) = \frac {4}{(a-b)^2} [-\frac{1}{6} (a-b)^3] = \frac{4}{6}(b-a) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.8)
 * }

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

$$ \int_{a}^{b}l_{2,2}(x) = \frac {2}{(b-a)^2} [-\frac{1}{12} (a-b)^3] = \frac{1}{6}(b-a) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.9)
 * }

With

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

$$ h=\frac {b-a}{2} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

Plug (Eq. 11.7)(Eq. 11.8)(Eq. 11.9) into (Eq. 11.3), yields

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

$$ I_2(f) = \int_{a}^{b} l_{0,2}(x)dx f(x_0) + \int_{a}^{b} l_{1,2}(x)dx f(x_1) + \int_{a}^{b} l_{2,2}(x)dx f(x_2) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 11.10)
 * }

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

$$ \Rightarrow I_2(f) = \frac {h}{3} (f(x_0) + 4f(x_1) + f(x_2)) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

Which is right the simple simpson's rule.

Problem
Given:

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

$$ f(x)=\frac {e^x-1}{x} on [-1,1], x_0=-1,x_n=1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.1)
 * }

Find: 1) Construct $$ f_n(x) $$ by lagrange interpolation functions for $$ n=1,2,4,8,16 $$ and then plot $$ f(x) $$ and $$ f_n(x) $$. 2) Compute weight functions for $$ n=1,2,4,8,16 $$ and then form tables. 3) Compute $$ I_n(x) $$ for $$ n=1,2,4,8 $$ and compare to $$ I $$ that is solved using WolframAlpha with more digits. 4) For $$ n=5 $$, plot $$ l_0,l_1,l_2 $$ and then predict the shapes of $$ l_3,l_4,l_5 $$ based on those.

Solution
1) According to the problem, we construct the fuctions from following equations:

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

$$ f_n(x)=p_n(x)=\sum_{i=0}^{n} l_{i,n}(x)f(x_i) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.2)
 * }
 * {| style="width:100%" border="0"

$$ l_{i,n}=\prod_{j=0,j \neq i}^{n} \frac{x-x_j}{x_i-x_j} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.3)
 * }

For n=1,

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

$$ x_0=-1,x_1=1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ l_{0,1}(x) = \frac {x-x_1}{x_0-x_1} = -\frac{1}{2}(x-1) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.4)
 * }

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

$$ l_{1,1}(x) = \frac {x-x_0}{x_1-x_0} = \frac{1}{2}(x+1) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.5)
 * }

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

$$ f(x_0) = \frac {e^{x_0}-1}{x_0} = 1-e^{-1} \displaystyle $$ $$ f(x_1) = \frac {e^{x_1}-1}{x_1} = e-1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ \Rightarrow f_1(x) = l_{0,1}(x)f(x_0) + l_{1,1}(x)f(x_1) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.6)
 * }

For n=2,

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

$$ x_0=-1,x_1=0,x_2=1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ l_{0,2}(x) = \frac {(x-x_1)(x-x_2)}{(x_0-x_1)(x_0-x_2)} = \frac{1}{2}x(x-1) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.7)
 * }

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

$$ l_{1,2}(x) = \frac {(x-x_0)(x-x_2)}{(x_1-x_0)(x_1-x_2)} = -(x+1)(x-1) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.8)
 * }

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

$$ l_{1,2}(x) = \frac {(x-x_0)(x-x_1)}{(x_2-x_0)(x_2-x_1)} = \frac{1}{2}(x+1)x $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.9)
 * }

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

$$ f(x_0) = \frac {e^{x_0}-1}{x_0} = 1-e^{-1} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ f(x_1) = \frac {e^{x_1}-1}{x_1} = 1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ f(x_2) = \frac {e^{x_2}-1}{x_2} = e-1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ \Rightarrow f_2(x) = l_{0,2}(x)f(x_0) + l_{1,2}(x)f(x_1) +  l_{2,2}(x)f(x_2)$$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.10)
 * }

For n=4,8,16, we use matlab for plotting with (Eq. 12.2) and (Eq. 12.3). The code is attached below:

The plots are shown as following:

2) For n=1, integrate previous (Eq. 12.4) and (Eq. 12.5) from [-1,1], we get,

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

$$ w_{0,1} = \int_{-1}^{1} l_{0,1}(x)dx = 1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ w_{1,1} = \int_{-1}^{1} l_{1,1}(x)dx = 1 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

For n=2, integrate previous (Eq. 12.7), (Eq. 12.8) and (Eq. 12.9) from [-1,1], we get,

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

$$ w_{0,2} = \int_{-1}^{1} l_{0,2}(x)dx = \frac{1}{3} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ w_{1,2} = \int_{-1}^{1} l_{1,2}(x)dx = \frac{4}{3} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ w_{2,2} = \int_{-1}^{1} l_{2,2}(x)dx = \frac{1}{3} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

For n=4,8,16, by apply LIET, we obtain the linear system equations:

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

$$ \sum_{i=0}^{n}w_{i,n}(x_i)^j = \frac{b^{j+1}-a^{j+1}}{j+1} $$ when $$ j=0,1,2,...,n $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.11)
 * }

Solving (Eq. 12.11) for n=4,8,16 by matlab and the code is attached below:

The results are shown in the following tables:

3) Integrate using Newton-Cotes Method, The integration could be expressed as:

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

$$ I_n(f) = \int_{a}^{b} f_n(x)dx = \sum_{i=0}^{n} w_{i,n}f(x_i) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.12)
 * }

With previously calculated value of weight functions, The results of the integrations are shown in the following table and compared with that obtained from WolframAlpha.

4) For n=5, the lagrange functions are computed from (Eq. 12.3) using matlab, the code is attached below:

The first three fuctions are plotted as follows:

Now, we will analyze the shapes of the lagrange functions. First it is obvious to identify:

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

$$ x_i = -x_{n-i} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

From (Eq. 12.3), we could derive

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

$$ l_{n-i,n}(x) = \prod_{j=0,j \neq n-i}^{n} \frac{x-x_j}{x_{n-i}-x_j} = \prod_{j=0,j \neq n-i}^{n} \frac{x+x_{n-j}}{x_{n-i}+x_{n-j}} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \prod_{j=0,n-j \neq i}^{n} \frac{-x-x_{n-j}}{x_{i}-x_{n-j}} = \prod_{k=0,k \neq i}^{n} \frac{-x-x_{k}}{x_{i}-x_{k}} = l_{i,n}(-x) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 12.13)
 * }

So from (Eq. 12.13), we can predict that $$ l_3 $$ and $$ l_2 $$, $$ l_4 $$ and $$ l_1 $$, $$ l_5 $$ and $$ l_0 $$ are axis-symmetric function pairs (With respect to y axis). We further plot $$ l_3 $$, $$ l_4 $$, $$ l_5 $$ in the following figure and the prediction is verified.

= 3.7 Integration Derivation for Error of Trap. rule =

Prove
Prove the following equation

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

$$ \frac{M_2}{2!} \int_{a}^{b} \left | (x-a)(x-b) \right | dx = \frac {(b-a)^3}{12} M_2 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.7.1)
 * }

Solution

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

$$ \frac{M_2}{2!} \int_{a}^{b} \left | (x-a)(x-b) \right | dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_2}{2!} \int_{a}^{b} (x-a)(b-x) dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_2}{2} \int_{a}^{b} [-x^2+(a+b)x-ab] dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_2}{2} [-\frac{1}{3}x^3+\frac{1}{2}(a+b)x^2-abx]\mid _{a}^{b}  $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_2}{2} [-\frac{1}{3}(b^3-a^3)+\frac{1}{2}(a+b)(b^2-a^2)-ab(b-a)] $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_2}{12} (b^3-3b^2a+3a^2b-a^3) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_2}{12} (b-a)^3 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

Prove complete. Egm6341.s11.team3.Xia 19:00, 11 February 2011 (UTC)

= 3.8 Integration Derivation for Error of Simple Simpson's rule =

Prove
Prove the following equation

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

$$ \frac{M_3}{3!} \int_{a}^{b} \left | (x-a)(x-\frac{a+b}{2})(x-b) \right | dx = \frac {(b-a)^4}{192} M_3 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.8.1)
 * }

Solution

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

$$ \frac{M_3}{3!} \int_{a}^{b} \left | (x-a)(x-\frac{a+b}{2})(x-b) \right | dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_3}{3!} [\int_{a}^{\frac{a+b}{2}} (x-a)(x-\frac{a+b}{2})(x-b) dx + \int_{\frac{a+b}{2}}^{b} (x-a)(x-\frac{a+b}{2})(b-x) dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_3}{3!} (\int_{a}^{\frac{a+b}{2}} - \int_{\frac{a+b}{2}}^{b}) [x^3-\frac{3}{2}(a+b)x^2+\frac{a^2+4ab+b^2}{2}x- \frac{a^2b+b^2a}{2}] dx $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_3}{3!} [\frac{1}{4}x^4-\frac{1}{2}(a+b)x^3+\frac{a^2+4ab+b^2}{4}x^2- \frac{a^2b+b^2a}{2}x] (\mid_{a}^{\frac{a+b}{2}} + \mid_{b}^{\frac{a+b}{2}}) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{2M_3}{6} [\frac{1}{64}(a+b)^4-\frac{1}{16}(a+b)(a+b)^3+\frac{a^2+4ab+b^2}{16}(a+b)^2- \frac{a^2b+b^2a}{4}(a+b)] $$ $$ - \frac{M_3}{6} [\frac{1}{4}a^4-\frac{1}{2}(a+b)a^3+\frac{a^2+4ab+b^2}{4}a^2- \frac{a^3b+b^2a^2}{2}] $$ $$ - \frac{M_3}{6} [\frac{1}{4}b^4-\frac{1}{2}(a+b)b^3+\frac{a^2+4ab+b^2}{4}b^2- \frac{a^2b^2+b^3a}{2}] $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_3}{192} (a^4-4a^3b-10a^2b^2-4ab^3+b^4) - \frac{M_3}{6} (-\frac{a^2b^2}{2}) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_3}{192} (a^4-4a^3b+6a^2b^2-4ab^3+b^4) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

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

$$ = \frac{M_3}{192} (b-a)^4 $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle
 * }

Prove complete. Egm6341.s11.team3.Xia 20:00, 11 February 2011 (UTC)

= 4.6 Fix Runge Phenomenon Use Non-uniform Distribution of Nodes =

Given:
Given a function

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

$$ f(x) = \frac{1}{1+4x^2}\quad $$ in interval $$ x \in [-5,5] \quad $$, $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 3.3.1)
 * }

Find:
1) Use the roots of legendre polynomials to compute $$\quad f_n^L(x) $$, plot $$\quad f_n^L(x) $$ and $$\quad f(x) $$, increase n till $$\quad I_n $$ is accurate to $$\quad O(10^{-6}) $$.

2) Use the roots of chebyshew polynomials to compute $$\quad f_n^L(x) $$, plot $$\quad f_n^L(x) $$ and $$\quad f(x) $$, Use Newton-Cotes to find $$\quad I_n $$ to $$\quad O(10^{-6}) $$ accuracy.

3) Use Gauss-legendre quadrature to find $$\quad I_n $$ to $$\quad O(10^{-6}) $$ accuracy.

4) Compare results of 1), 2) and 3).

Solution:
For comparison, the integration result is computed from WolframAlpha:

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

$$ \int_{-5}^{5} \frac{1}{1+4x^2} dx = 1.471127674303734 \quad $$, $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.1)
 * }

1) To apply legendre polynomials, we convert the variable domain from [-5 5] to [-1 1] by the following transformation:

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

$$ f^*(X) = \frac{1}{1+100X^2}\quad $$ and $$ X=x/5 \in [-1,1] \quad $$, $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.1)
 * }

And the integration realtions are given by:

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

$$ I(x) = \int_{-5}^{5} f(x)dx = \int_{-1}^{1} f^*(X)d(5X) = 5I^*(X) $$, $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.2)
 * }

The legendre polynomial are given by:
 * {| style="width:100%" border="0"

$$ P_{n}(X) = \sum_{i=0}^{[n/2]} (-1)^i\frac{(2n-2i)!X^{n-2i}}{2^ni!(n-i)!(n-2i)!} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.3)
 * }

From Newton-Cotes method:

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

$$ f_n^{*L}(X)=\sum_{i=0}^{n+1} l_{i,n+1}(X)f^*(X_i) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.4)
 * }
 * {| style="width:100%" border="0"

$$ l_{i,n+1}=\prod_{j=0,j \neq i}^{n+1} \frac{X-X_j}{X_i-X_j} $$ $$ where $$ {X_i, i=1,2,...,n} \quad $$ are the roots of (Eq. 4.6.1), $$ X_0=-1, X_{n+1}=1\quad $$. $$ I_n(X) \quad $$ is given by:
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.5)
 * }

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

$$ I_{n+1}^*(X)=\sum_{i=0}^{n+1} w_{i,n+1} f^*(X_i) $$ $$ the $$ w_{i,n+1} \quad $$ is estimated by solving the linear system matrix:
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.6)
 * }

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

$$ \sum_{i=0}^{n+1} w_{i,n+1}(X_i)^j = \frac {b^{j+1}-a^{j+1}}{j+1} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.7)
 * }

With matlab solving the roots for (Eq. 4.6.1), we can obtain $$\quad f_n^L(x) $$ and $$\quad I_n $$. The code is attached below:

The resultant plots are shown below:

The Integrations and Errors are listed in the table below:

2) To apply chebyshew polynomials, we convert the variable domain from [-5 5] to [-1 1] by the same transformation as (Eq. 4.6.1) and (Eq. 4.6.2):

The roots of chebyshew polynomials are given by:
 * {| style="width:100%" border="0"

$$ X_{i,n} = cos(i\frac{\pi}{n})$$ and $$ i=0,1,2,...,n \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.8)
 * }

Then similar to the previous section, we apply Newton-Cotes method and lagrange interpolation method as in (Eq. 4.6.4) and (Eq. 4.6.6). Using matlab to obtain $$\quad f_n^L(x) $$ and $$\quad I_n $$. The code is attached below:

The resultant plots are shown below:

The Integrations and Errors are listed in the table below:

3) To apply Gauss-legendre quadrature, we convert the variable domain from [-5 5] to [-1 1] by the same transformation as (Eq. 4.6.1) and (Eq. 4.6.2).

The legendre polynomial are given by (Eq. 4.6.3), from which the roots are obtained. The weight functions are computed by:

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

$$ w_i= \frac {-2}{(n+1)P_n^{(1)}(X_i)P_{n+1}(X_i)} $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.9)
 * }

The integration $$\quad I_n^* $$ is therefore given by:

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

$$ I_n^* = \sum_{i=1}^{n} w_i f^*(X_i) $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 4.6.10)
 * }

Using matlab to solve for the integrations. The code is attached below:

The Integrations and Errors are listed in the table below:

4) Explanation for not achieving $$ O(10^{-6}) \quad $$: Basically, in this assignment, n is computed to 32. For higher n(such as 64), the matlab "roots" function would fail and generated complex values for legendre poly roots. Moreover, the matrix for computing the weights also became "ill" so that it is considered "single" by the machine. Therefore, n=32 is the best accuracy we got for this problem.

Comparing the results of the error, it is obvious that Newton-Cotes method (first 2 cases) has a better accuracy than Gauss-Legendre Quadrature method when n is small, while the Gauss-Legendre Quadrature method converges faster than the Newton-Cotes method. Comparing the different interpolation approaches for the Newton-Cotes method, they have similar behaviors of accuray. When n=32, all the three methods are of the order of $$ O(10^{-3}) \quad $$.

Egm6341.s11.team3.Xia 5:00, 1 March 2011 (UTC)

Results:
The integration results for a), c) and d) are shown and compared in the following table:

The Romberg table for b) is:

(2) Circumference for Ellipse: Method1
From definition of eccentricity of ellipse :

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

$$ e = \sqrt {1-\frac{b^2}{a^2}} = \frac{\sqrt{91}}{10} \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 6.7.11)
 * }

From (Eq. 6.7.2), (Eq. 6.7.3), we derive the Integration as:

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

$$ C = I = \int_{0}^{2 \pi} a(1-e^2)(1-ecos \theta)^{-2}\sqrt{(1 +e^2 - 2e cos \theta)} d \theta \quad $$ $$ = \int_{0}^{2 \pi} \frac{9}{10} (1- \frac{\sqrt{91}}{10} cos \theta)^{-2}(\frac{191}{100} - \frac{\sqrt{91}}{5} cos \theta)^{\frac{1}{2}} d \theta \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 6.7.12)
 * }

Using WolframAlpha, the exact solution is calculated to be :

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

$$ C_{exact} = I_{exact} = 43.8591 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 6.7.13)
 * }

Following the same approaches as in (Eq. 6.7.7), (Eq. 6.7.8) and (Eq. 6.7.9), the integration is computed using the following code:

Results:
The integration results for a), c) and d) are shown and compared in the following table:

The Romberg table for b) is:

(3) Circumference for Ellipse: Method2
The 2nd method of circumference for ellipse follows (Eq. 6.7.4). Similar to previous calculations, the function integrated here is defined as:

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

$$ f(\theta) = 4a[1-e^2 sin^2 \theta]^{\frac{1}{2}} = 40\sqrt{1- 0.91 sin^2 \theta} \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 6.7.14)
 * }

Using WolframAlpha, the exact solution is calculated to be :

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

$$ C_{exact} = I_{exact} = 43.8591 \quad $$ $$
 * style="width:95%" |
 * style="width:95%" |
 *
 $$ \displaystyle (Eq. 6.7.13)
 * }

Which is identical to (Eq. 6.7.13), following the same approaches as in (Eq. 6.7.7), (Eq. 6.7.8) and (Eq. 6.7.9), the integration is computed using the following code:

Results:
The integration results for a), c) and d) are shown and compared in the following table:

The Romberg table for b) is:

Egm6341.s11.team3.Xia 18:00, 31 March 2011 (UTC)

=References=

title: User:EML4500.f08.JAMAMA/NM/Xia


text:
User:EML4500.f08.JAMAMA/NM/elango

= 3.3.1 Plotting figures similar to proof of LIET =

Given
Given a function

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

$$ f(x) = \frac{exp(x)-1}{x} \quad $$ in interval $$ x \in [-1, 1] \quad $$, $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 3.3.1)
 * }

Find
Find by Lagrange polynomial,weights and $$ q(n+1) \quad $$ for, $$ L_2,4^L(x) \quad $$, at 5 nodes. Nodes were chosen to be equidistant $$ x_0 = -1, x_1 = \frac{-1}{2}, x_2 = \frac 0, x_3 = \frac{1}{2}, x_4 = 1 $$.

also find errors at points $$ x = 0.5 \quad $$ and $$ t = 5 \quad $$

Solution
Consider a Lagrange method on approximating exact function: $$f_n^L\left( x \right)=\,\sum_{i=0}^{n}l_{i,n}\left ( x \right)f\left ( x_i \right)$$

For given 5 nodes at $$ n = 4 \rightarrow x_{0} = -1, x_{1} = -0.5, x_{2} = 0, x_{3} = 0.5, x_{4} = 1 $$.

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

f_{4}^L=\sum_{i=0}^{4}l_{i,4}(x)f(x_{i})=l_{0,4}(x)f(x_{0})+l_{1,4}(x)f(x_{1})+l_{2,4}(x)f(x_{2})+l_{3,4}(x)f(x_{3})+l_{4,4}(x)f(x_{4}) $$ $$ where
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle (Eq. 3.3.1.3)
 * }
 * }
 * {| style="width:100%" border="0" align="left"

l_{i,n}=\prod_{j=0,j\neq 0}^{4}\frac{x-x_{j}}{x_{i}-x_{j}} $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * {| style="width:100%" border="0" align="left"
 * {| style="width:100%" border="0" align="left"

l_{0,4}=\prod_{j=0,j\neq 0}^{4}\frac{x-x_{j}}{x_{0}-x_{j}} $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * {| style="width:100%" border="0" align="left"
 * {| style="width:100%" border="0" align="left"

l_{1,4}=\prod_{j=0,j\neq 1}^{4}\frac{x-x_{j}}{x_{1}-x_{j}} $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * {| style="width:100%" border="0" align="left"
 * {| style="width:100%" border="0" align="left"

l_{2,4}=\prod_{j=0,j\neq 2}^{4}\frac{x-x_{j}}{x_{2}-x_{j}} $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * {| style="width:100%" border="0" align="left"
 * {| style="width:100%" border="0" align="left"

l_{3,4}=\prod_{j=0,j\neq 3}^{4}\frac{x-x_{j}}{x_{3}-x_{j}} $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * {| style="width:100%" border="0" align="left"
 * {| style="width:100%" border="0" align="left"

l_{4,4}=\prod_{j=0,j\neq 4}^{4}\frac{x-x_{j}}{x_{4}-x_{j}} $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

For example:

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

l_{0,4}=\prod_{j=0,j\neq 0}^{4}\frac{x-x_{j}}{x_{0}-x_{j}} = \frac{x-x_1}{x_{0}-x_{1}}\frac{x-x_{2}}{x_{0}-x_{2}}\frac{x-x_{3}}{x_{0}-x_{3}}\frac{x-x_{4}}{x_{0}-x_{4}} = \left(\frac{x-(-0.5)}{(-1)-(-0.5)} \right) \left(\frac{x-0}{(-1)-0} \right ) \left( \frac{x-0.5}{(-1)-0.5} \right) \left (\frac{x-1}{(-1)-1} \right ) $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Therefore, plotting $$ Eq. (3.3.1.3) \quad $$ with given parameters we get:

Maple code was written to plot the graphs and is given below, and the equations used are provided.

q_{n+1}(x) := (x - x_0)(x - x_1)(x - x_2)(x - x_3)\cdots(x - x_n) $$ $$
 * $$\displaystyle
 * $$\displaystyle (Eq. 3.3.1.2)

Maple code was written to plot the graphs and is given below, and the equations used are provided.

$$ |e_4^L(t)| = 0.0046570,|e_4^L(x)| = 0.3955996 \quad $$

= 3.3.2 Runge-Phenomenon =

Given
Given a function

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

$$ f(x) = \frac{1}{1+4x^2}\quad $$ in interval $$ x \in [-5,5] \quad $$, $$
 * style="width:95%" |
 * style="width:95%" |
 *  $$ \displaystyle (Eq. 3.3.2)
 * }

Find
Find the error estimate for the function given above and plot, $$ L_3,8^L(x) \quad $$, for n=8. Nodes were chosen to be equidistant $$ a = -5,b = 5 \quad$$.

Solution
When n = 8 $$ \rightarrow x_{0} = -5, x_{1} = -3.75, x_{2} = -2.5, x_{3} = -1.25, x_{4} = 0, x_{5} = 1.25, x_{6} = 2.5, x_{7} = 3.75, x_{8} = 5 $$.

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

 * $$\displaystyle
 * $$\displaystyle

f_{8}=\sum_{i=0}^{8}l_{i,8}(x)f(x_{i})=l_{0,8}(x)f(x_{0})+l_{1,8}(x)f(x_{1})+l_{2,8}(x)f(x_{2})+l_{3,8}(x)f(x_{3})+l_{4,8}(x)f(x_{4})+l_{5,8}(x)f(x_{5})+l_{6,8}(x)f(x_{6})+l_{7,8}(x)f(x_{7})+l_{8,8}(x)f(x_{8}) $$ $$
 * $$\displaystyle (Eq.3.3.2.1)
 * }
 * }

where,

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

l_{3,8}=\prod_{j=0,j\neq 3}^{8}\frac{x-x_{j}}{x_{3}-x_{j}} $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Maple code was written to solve for the error estimates and to plot the req function

$$ |e_8^L(t)| = 0.001415387826,|e_8^L(x)| = -0.2033016134e-3 \quad $$

Ref

 * }

Given that, $$ i=1 \,$$

$$ \frac{c_{1}}{3!}+c_{3}=0\, $$

$$ p_{1}(t)=-t,i.e.,c_{1}=-1 \, $$

$$ \Rightarrow c_{3}=\frac{1}{6}\, $$

when $$ i=2 \,$$

$$ \frac{c_{1}}{5!}+\frac{c_{3}}{3!}+c_{5}=0\, $$

$$ c_{3}=\frac{1}{6}\, $$

$$ \Rightarrow c_{5}=-\frac{7}{360}\, $$

when $$ i=3 \,$$

$$ \frac{c_{1}}{7!}+\frac{c_{3}}{5!}+\frac{c_{5}}{3!}+c_{7}=0\, $$

$$ c_{5}=-\frac{7}{360}\, $$

$$ \Rightarrow c_{7}=\frac{31}{15120}\, $$

when $$ i=4 \,$$

$$ \frac{c_{1}}{9!}+\frac{c_{3}}{7!}+\frac{c_{5}}{5!}+\frac{c_{7}}{3!}+c_{9}=0\, $$

substitute the values and solve for c_{9} ,

$$ \Rightarrow c_{9}=-\frac{127}{604800}\, $$

Now we calculate the pair polynomials using the expressions given below for odd and even terms, and then the constants are substituted in them ,

$$ p_{2i}(t)=\sum_{j=0}^{i}c_{2j+1}\frac{t^{2(i-j)}}{[2(i-j)]!}\, $$

$$ p_{2i+1}(t)=\sum_{j=0}^{i}c_{2j+1}\frac{t^{2(i-j)+1}}{[2(i-j)+1]!}+c_{2i+2}\, $$

The expressions are generated as follows,

Now we get the final Pair Polynomials ,by just substituting the constants in the maple generated equations and the pairs are listed below :

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

 * style="width:50%; padding:10px; border:2px solid #8888aa" |
 * $$ \displaystyle
 * $$ \displaystyle

p_{2}=-\frac{t^{2}}{2}+\frac{1}{6},

p_{3}=-\frac{t^{3}}{6}+\frac{1t}{6}\,$$ $$
 * $$\displaystyle
 * }
 * }

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

 * style="width:50%; padding:10px; border:2px solid #8888aa" |
 * $$ \displaystyle
 * $$ \displaystyle

p_{4}=-\frac{t^{4}}{24}+\frac{t^{2}}{12}-\frac{7}{360},

p_{5}=-\frac{t^{5}}{120}+\frac{t^{3}}{36}-\frac{7t}{360}\,$$ $$
 * $$\displaystyle
 * }
 * }

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

 * style="width:50%; padding:10px; border:2px solid #8888aa" |
 * $$ \displaystyle
 * $$ \displaystyle

p_{6}=-\frac{t^{6}}{720}+\frac{t^{4}}{144}-\frac{7t^{2}}{720}+\frac{31}{15120},

p_{7}=-\frac{t^{7}}{5040}+\frac{t^{5}}{720}-\frac{7t^{3}}{2160}+\frac{31t}{15120}\,$$ $$
 * $$\displaystyle
 * }
 * }

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

 * style="width:50%; padding:10px; border:2px solid #8888aa" |
 * $$ \displaystyle
 * $$ \displaystyle

p_{8}=-\frac{t^{8}}{40320}+\frac{t^{6}}{4320}-\frac{7t^{4}}{8640}+\frac{31t^{2}}{30240}-\frac{127}{604800},

p_{9}=-\frac{t^{9}}{362880}+\frac{t^{7}}{30240}-\frac{7t^{5}}{43200}+\frac{31t^{3}}{90720}-\frac{127t}{604800}\,$$ $$
 * $$\displaystyle
 * }
 * }

title: User:EML4500.f08.JAMAMA/NM/elango


text:
User:EML4500.f08.JAMAMA/NM/JA

= 7.1 =

Given: Verhulst Equation with initial conditions
Verhulst equation:

$$ \frac{d}{dt}x(t) = rx(1 - \frac{x}{x_{max}}) \ $$

$$ x_{max} = 15 \ $$ and $$ r = 1.4 \ $$

with initial condition of $$ x_0 \ $$

for time $$ t \in [0, 10] \ $$

Case 1:

$$ x_0 = 3 < \frac{1}{2} x_{max} \ $$

Case 2:

$$ x_0 = 9 > \frac{1}{2} x_{max} \ $$

Find: By using Hermite - Simpson rule, integrate given Verhulst equation
From given Prof. Vu Quoc's notes for

$$ z_{i+ 1} = z_i + \frac{h/2}{3} \left ( f_i + 4f_{i+\frac{1}{2}}[g(z_i,z_{i+1})],t_{i+\frac{1}{2}} + f_{i+1} \right ) \ $$

Solution
From given Prof. Vu Quoc's notes for

$$ z_{i+ 1} = z_i + \frac{h/2}{3} \left ( f_i + 4f_{i+\frac{1}{2}}[g(z_i,z_{i+1})],t_{i+\frac{1}{2}} + f_{i+1} \right ) \ $$

,where

$$ f_i = f(z_i, f_i) \ $$

$$ f_{i+1} = f(z_{i+1}, t_{i+1}) \ $$

$$ f_{i+\frac{1}{2}} = f[g(z_i,z_{i+1},t_{i+\frac{1}{2} })] = f(z_{i+\frac{1}{2}}, t_{i+\frac{1}{2}}) \ $$

,where $$ t_{i+\frac{1}{2}} = t_i + \frac{h}{2} \ $$ and as it was proven in HW 6.6 as well as given by Prof. Vu Quoc's

$$ z_{i+ \frac{1}{2}} = \frac{1}{2}(z_i + z_{i+1}) + \frac{h}{8}(f_i - f_{i+1}) \ $$

,note that $$ h \ $$ is the time step

Also specifying the equation to be of a form $$ F(z_{i+ 1})= 0 \ $$

$$ F(z_{i+ 1})= -z_{i+ 1} + z_i + \frac{h/2}{3} \left ( f_i + 4f([g(z_i,z_{i+1})],t_{i+\frac{1}{2}}) + f_{i+1} \right ) = 0\ $$

and, finally, as mentioned in

$$ z_{i+1}^{(k+1)} = z_{i+1}^{(k)} - \left [ \frac{d}{dz}F(z_{i+1}^{(k)}) \right ]^{-1} F(z_{i+1}^{(k)}) \ $$

in our case, for a given function, we have first order derivative. Therefore

$$ z_{i+1}^{(1)} = z_{i+1} - \left [ \frac{d}{dz}F(z_{i+1}) \right ]^{-1} F(z_{i+1}) \ $$

Then check the condition for convergence as it is shown in

$$ z_{i+1}^{(1)} - z_{i+1} ?< Tolerance \ $$

If condition is not met, choose $$ z_{i+1}^{(1)} = z_{i+2} \ $$ and iterate again to get $$ z_{i+2}^{(1)} \ $$ until the tolerance is met

Egm6341.s11.team3.rakesh 02:24, 21 April 2011 (UTC)

title: User:EML4500.f08.JAMAMA/NM/JA


text:
Template:MyProjectNav

This template was transcluded from Wikieducator - Template:MyProjectNav

title: Template:MyProjectNav


text:
User:Egm6341.s11.team2.zou/Team Format

=Problem 2.1: "CONCISELY STATE PROBLEM OBJECTIVE HERE--CONCISELY"= NOTHING SHOULD GO HERE

Given: OR "PROBLEM STATEMENT:" OR "BACKGROUND" OR "THEOREM"
THIS SHOULD INCLUDE EITHER WHAT IS EXPLICITLY GIVEN IN CLASS, BACKGROUND INFORMATION SUCH AS THEOREMS OR LEMMAS, AND A FIGURE OR GRAPH IF IT PERTAINS TO THE PROBLEM STATEMENT

IF EQUATIONS ARE USED HERE AND TO BE REFERENCE LATER IN THE SOLUTION PLEASE USE CENTER JUSTIFIED EQUATIONS WITH THE FOLLOWING LABELING

Find: OR "PROVE" OR "DETERMINE" OR "EVALUATE"
IF THIS INCLUDES SUBPARTS LET TRY TO USE THE FOLLOWING CONVENTION A) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$ B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$ C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$ ETC...

"SOLUTION" OR "PROOF"
A) COPY THE STATEMENT "DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$" SO THE TA KNOWS WHICH PART HE IS GRADING WORDS! USE LOTS OF WORD TO DESCRIBE THE MATHEMATICS

CONCLUDE WITH THEREFORE AND BOX SOLUTION

B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$

C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$

=Problem 2.1: "CONCISELY STATE PROBLEM OBJECTIVE HERE--CONCISELY"= NOTHING SHOULD GO HERE

Given: OR "PROBLEM STATEMENT:" OR "BACKGROUND" OR "THEOREM"
THIS SHOULD INCLUDE EITHER WHAT IS EXPLICITLY GIVEN IN CLASS, BACKGROUND INFORMATION SUCH AS THEOREMS OR LEMMAS, AND A FIGURE OR GRAPH IF IT PERTAINS TO THE PROBLEM STATEMENT

IF EQUATIONS ARE USED HERE AND TO BE REFERENCE LATER IN THE SOLUTION PLEASE USE CENTER JUSTIFIED EQUATIONS WITH THE FOLLOWING LABELING

Find: OR "PROVE" OR "DETERMINE" OR "EVALUATE"
IF THIS INCLUDES SUBPARTS LET TRY TO USE THE FOLLOWING CONVENTION A) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$ B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$ C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$ ETC...

"SOLUTION" OR "PROOF"
A) COPY THE STATEMENT "DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$" SO THE TA KNOWS WHICH PART HE IS GRADING WORDS! USE LOTS OF WORD TO DESCRIBE THE MATHEMATICS

CONCLUDE WITH THEREFORE AND BOX SOLUTION

B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$

C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$

Given: OR "PROBLEM STATEMENT:" OR "BACKGROUND" OR "THEOREM"
THIS SHOULD INCLUDE EITHER WHAT IS EXPLICITLY GIVEN IN CLASS, BACKGROUND INFORMATION SUCH AS THEOREMS OR LEMMAS, AND A FIGURE OR GRAPH IF IT PERTAINS TO THE PROBLEM STATEMENT

IF EQUATIONS ARE USED HERE AND TO BE REFERENCE LATER IN THE SOLUTION PLEASE USE CENTER JUSTIFIED EQUATIONS WITH THE FOLLOWING LABELING

Find: OR "PROVE" OR "DETERMINE" OR "EVALUATE"
IF THIS INCLUDES SUBPARTS LET TRY TO USE THE FOLLOWING CONVENTION A) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$ B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$ C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$ ETC...

"SOLUTION" OR "PROOF"
A) COPY THE STATEMENT "DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$" SO THE TA KNOWS WHICH PART HE IS GRADING WORDS! USE LOTS OF WORD TO DESCRIBE THE MATHEMATICS

CONCLUDE WITH THEREFORE AND BOX SOLUTION

B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$

C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$

Given: OR "PROBLEM STATEMENT:" OR "BACKGROUND" OR "THEOREM"
THIS SHOULD INCLUDE EITHER WHAT IS EXPLICITLY GIVEN IN CLASS, BACKGROUND INFORMATION SUCH AS THEOREMS OR LEMMAS, AND A FIGURE OR GRAPH IF IT PERTAINS TO THE PROBLEM STATEMENT

IF EQUATIONS ARE USED HERE AND TO BE REFERENCE LATER IN THE SOLUTION PLEASE USE CENTER JUSTIFIED EQUATIONS WITH THE FOLLOWING LABELING

Find: OR "PROVE" OR "DETERMINE" OR "EVALUATE"
IF THIS INCLUDES SUBPARTS LET TRY TO USE THE FOLLOWING CONVENTION A) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$ B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$ C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$ ETC...

"SOLUTION" OR "PROOF"
A) COPY THE STATEMENT "DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 1 $$" SO THE TA KNOWS WHICH PART HE IS GRADING WORDS! USE LOTS OF WORD TO DESCRIBE THE MATHEMATICS

CONCLUDE WITH THEREFORE AND BOX SOLUTION

B) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 2 $$

C) DETERMINE $$\displaystyle p_{n} $$ SUCH THAT $$\displaystyle n = 3 $$

=Contributing Members=

title: User:Egm6341.s11.team2.zou/Team Format


text:
User:Oannasoro

ro

title: User:Oannasoro


text:
Template:Vbar

|

title: Template:Vbar


text:
Category:Resources by who participants are

This category should be used to guide participants to the help/introductory resources they may need depending on why they have chosen to participate at Wikiversity. Resources should be added to subcategories and not directly to this category.


text:
SPIR608/Clubs

This is a list of some of the game clubs in London

Tuesdays

 * Bush Bash, Shepherds Bush or Gold Hawk Road, 7pm

Wednesdays

 * London On Board, St Paul's, 5:45pm
 * Swiggers, London Bridge, 6pm

Sundays
Chestnut Lodge Wargames group
 * Sunday 6th February: Balkans Political Game and Goodwood at Anerley Methodist Church
 * Sunday 6th March: Hell of the North and Carry on Colonelsat Anerley Methodist Church
 * Sunday 3rd April: (Continued from Saturday 2nd April) Strategic Triad, Seven Years War Strategic, Once More over the Hill and 1814 Campaign

title: SPIR608/Clubs


text:
SPIR608/Weblinks

General Websites

 * BoardGameGeek: reviews, photographs and forums for lovers of board games.
 * Brettspiel
 * Journal of Board Game Design
 * Mike Doyle, The Function(s) of Game Art
 * Computer Games Design at UEL

Scott Nicholson
Scott Nicholson is a an Associate Professor at Syracuse University
 * Board Games with Scott: Video reviews of loads of board games.
 * GAMBIT Class: CMS.608 Scott Nicholson on "Board Game Publishing"

title: SPIR608/Weblinks


text:
SPIR608/Components

On-line

 * Spiel Material
 * Meeple
 * Kallistra
 * War Bases
 * Board Game Extras

title: SPIR608/Components


text: Category:Exploring Wikiversity as a community builder