User:Egm6341.s10.Team4.yunseok/HW22

= Problem 5: Comparison of Error between Taylor series and Newton-Cotes formula =

Given
Refer Lecture slide 11-1 for problem statement
 * {| style="width:100%" border="0" align="left"

f(x) = \sin(x) $$ $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle x \in [0,\pi] $$
 * $$\displaystyle \longrightarrow (1)


 * }.
 * }.

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

\left | f^T_n\Big(\frac{7\pi}{8}\Big)-f\Big(\frac{7\pi}{8}\Big) \right | \le \underbrace{\left | f^L_4\Big(\frac{7\pi}{8}\Big)-f\Big(\frac{7\pi}{8}\Big) \right |}_{\left | E^L_4(\frac{7\pi}{8})\right |} \le \frac{\left | q_5(\frac{7\pi}{8}) \right |}{5!} $$ $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle \longrightarrow (2)


 * }.
 * }.

Solution
$$\displaystyle f^L_4\Big(\frac{7\pi}{8}\Big)$$, $$\displaystyle \left | E^L_4(\frac{7\pi}{8})\right |$$, $$\displaystyle \frac{\left | q_5(\frac{7\pi}{8}) \right |}{5!}$$, is claculated using Matlab.


 * $$\displaystyle f^L_4\Big(\frac{7\pi}{8}\Big) = \sum_{i=0}^4 l_{i,4}\Big(\frac{7\pi}{8}\Big)\sin(x_i)= 0.3812$$
 * $$\displaystyle \left | E^L_4(\frac{7\pi}{8})\right | = 0.0014808$$
 * $$\displaystyle \frac{\left | q_5(\frac{7\pi}{8}) \right |}{5!} = 0.0081716$$

At first, we can confirm that the seconde inequality of Equation(2)


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

\left | E^L_4(\frac{7\pi}{8})\right | = 0.0014808 \le \displaystyle \frac{\left | q_5(\frac{7\pi}{8}) \right |}{5!} = 0.0081716 $$ $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle \longrightarrow (3)


 * }.
 * }.

Then, the Taylor series were calculated until the Error or taylor become smaller than that of Newton-Cotes formula at n=4. the result is like that

$$ \begin{array}{|c||c|c|} n & E^T_n(t) & E^T_n(t) - E^L_4(t) \\ \hline 1&1.7128&1.7113\\ 2&0.34976&0.3483\\	3&0.54236&0.5409\\ 4&0.10444&0.1030\\ 5&0.067533&0.0661\\	6&0.011256&0.0098\\	7&0.0045295&0.0030\\	8&0.00065513&-0.0008\\ \hline \end{array} $$

Matlab Code:

= Problem 10: Proof of Error =

Given
Refer Lecture slide 13-2 for problem statement
 * {| style="width:100%" border="0" align="left"

\begin{align} n=2, \\ & q_3(x)=(x-x_o)(x-x_1)(x-x_2)\\ & x_o=a, x_2=b, x_1=\frac{a+b}{2} \end{align} $$
 * $$\displaystyle
 * $$\displaystyle


 * }.
 * }.

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

\left | E_2 \right | = \frac{M_3}{3!}\int_{a}^{b}\left|(x-x_o)(x-x_1)(x-x_2) \right |\, dx = \frac{(b-a)^4}{192}M_3 $$ $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle \longrightarrow (1)
 * }.
 * }.

Solution

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

\begin{align} \int_{a}^{b}\left|(x-a)(x-\frac{a+b}{2})(x-b) \right |\, dx &= 2\int_{a}^{\frac{a+b}{2}}(x-a)\Big(x-\frac{a+b}{2}\Big)(x-b)\,dx\\ &=2\int_{a}^{\frac{a+b}{2}}x^3-\frac{3(a+b)}{2}x^2+\frac{1}{2}(a^2+4ab+b^2)x - \frac{1}{2}ab(a+b)\, dx\\ &=\frac{1}{2}x^4-(a+b)x^3+\frac{1}{2}(a^2+4ab+b^2)x^2 - ab(a+b)x\, \Big |^{\frac{a+b}{2}}_a\\ &=\underbrace{\frac{1}{2}\Bigg(\Big(\frac{a+b}{2}\Big)^4-a^4\Bigg)}_{1st}-\underbrace{(a+b)\Bigg(\Big(\frac{a+b}{2}\Big)^3-a^3\Bigg)}_{2nd}+\underbrace{\frac{1}{2}(a^2+4ab+b^2)\Bigg(\Big(\frac{a+b}{2}\Big)^2-a^2\Bigg)}_{3rd} - \underbrace{ab(a+b)\Big(\frac{a+b}{2}-a\Big)}_{4th} \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * }.

1st and 3rd term can be bundled by $$\displaystyle \frac{1}{2}\Bigg(\Big(\frac{a+b}{2}\Big)^2-a^2\Bigg)$$


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

\begin{align} &=\frac{1}{2}\Bigg(\Big(\frac{a+b}{2}\Big)^2-a^2\Bigg)\Bigg(\Big(\frac{a+b}{2}\Big)^2+a^2+(a^2+4ab+b^2)\Bigg)\\ &=\frac{1}{2}\Bigg(\Big(\frac{a+b}{2}\Big)+a\Bigg)\Bigg(\Big(\frac{a+b}{2}\Big)-a\Bigg)\Big(\frac{a^2+2ab+b^2}{4}+2a^2+4ab+b^2\Bigg)\\ &=\frac{1}{32}(3a+b)(b-a)(9a^2+18ab+5b^2)\\ \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * {| style="width:１０0%" border="0" align="left"
 * {| style="width:１０0%" border="0" align="left"

$$\displaystyle =\frac{1}{32}(b-a)(27a^3+63a^2b+33ab^2+5b^3) $$ $$
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * $$\displaystyle \longrightarrow (2)
 * style= |
 * }.
 * }.

2nd and 4th term can be bundled by $$\displaystyle -(a+b)\Big(\frac{a+b}{2}-a\Big)$$
 * {| style="width:100%" border="0" align="left"

\begin{align} &=-(a+b)\Big(\frac{a+b}{2}-a\Big)\Bigg(\Big(\frac{a+b}{2}\Big)^2+\frac{a+b}{2}a + a^2 + ab \Bigg)\\ &=-\frac{1}{8}(a+b)(b-a)(7a^2+8ab+b^2)\\ &=-\frac{1}{8}(b-a)(7a^3+15a^2b+9ab^2+b^3)\\ \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * {| style="width:100%" border="0" align="left"

$$\displaystyle =-\frac{1}{32}(b-a)(28a^3+60a^2b+36ab^2+4b^3) $$ $$
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * $$\displaystyle \longrightarrow (3)
 * style= |
 * }.
 * }.

Then, (2) + (3),
 * {| style="width:100%" border="0" align="left"

\begin{align} &\int_{a}^{b}\left|(x-a)(x-\frac{a+b}{2})(x-b) \right |\, dx\\ &=\frac{1}{32}(b-a)(27a^3+63a^2b+33ab^2+5b^3)-\frac{1}{32}(b-a)(28a^3+60a^2b+36ab^2+4b^3)\\ &=\frac{1}{32}(b-a)\Big((27a^3+63a^2b+33ab^2+5b^3)-(28a^3+60a^2b+36ab^2+4b^3)\Big)\\ &=\frac{1}{32}(b-a)(-a^3+3a^2b-3ab^2+b^3)\\ &=\frac{1}{32}(b-a)(b-a)^3\\ \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * {| style="width:100%" border="0" align="left"
 * {| style="width:100%" border="0" align="left"

$$\displaystyle = \frac{1}{32}(b-a)^4 $$ $$
 * style="width:5%; padding:10px; border:2px solid #8888aa" |
 * style="width:5%; padding:10px; border:2px solid #8888aa" |
 * $$\displaystyle \longrightarrow (4)
 * style= |
 * }.
 * }.

So, plug in (4) to (1) equation.


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

$$\displaystyle = \left | E_2 \right | = \frac{M_3}{3!}\int_{a}^{b}\left|(x-x_o)(x-x_1)(x-x_2) \right |\, dx = \frac{M_3}{3!}\Bigg(\frac{1}{32}(b-a)^4\Bigg) = \frac{(b-a)^4}{192}M_3 $$
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * style= |
 * }.
 * }.

= Problem 11: Comparison of I and $$I_2$$of 3rd polynomial equation =

Given
Refer Lecture slide 13-3 for problem statement
 * {| style="width:70%" border="0" align="left"

f(x)=6x^3-2x^2+8x+3 $$ x \in [0,1] $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle


 * }.
 * }.

Find
Find $$\displaystyle I, I_2$$ (simpson's rule) and compare them.


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

\begin{align} &I=\int_{0}^{1} 6x^3-2x^2+8x+3,\ dx\\ &I_2=\frac{h}{3}\Big[f(x_0)+4f(x_1)+f(x_2)\Big] \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * }.

Solution
$$\displaystyle I $$ is integrated like that
 * {| style="width:100%" border="0" align="left"

\begin{align} I&=\int_{0}^{1} 6x^3-2x^2+8x+3,\ dx\\ &= \frac{6}{4}x^4-\frac{2}{3}x^3+\frac{8}{2}x^2+3x,\ Big|^1_0 \\ &= \frac{3}{2}(1^4-0^4)-\frac{2}{3}(1^3-0^3)+4(1^2-0^2)+3(1-0)\\ &= \frac{3}{2}(1^4-0^4)-\frac{2}{3}(1^3-0^3)+4(1^2-0^2)+3(1-0)\\ \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * {| style="width:100%" border="0" align="left"
 * {| style="width:100%" border="0" align="left"

$$\displaystyle I= \frac{3}{2}-\frac{2}{3}+4+3= \frac{47}{6} $$ $$
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * $$\displaystyle \longrightarrow (1)
 * style= |
 * }.
 * }.

$$\displaystyle I_2 $$ is,
 * {| style="width:100%" border="0" align="left"

I_2=\frac{h}{3}\Big[f(x_0)+4f(x_1)+f(x_2)\Big] $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * }.

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

\begin{align} &h=\frac{1}{2}\\ &f(x_0)=f(0)=3\\ &f(x_1)=f(\frac{1}{2})=6(\frac{1}{2})^3-2(\frac{1}{2})^2+8(\frac{1}{2})+3 = \frac{6}{8} - \frac{2}{4} + \frac{8}{2}+3 = \frac{29}{4}\\ &f(x_2)=f(1)=6(1)^3-2(1)^2+8(1)+3 =6-2+8+3 = 15\\ \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * }.

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

$$\displaystyle I_2=\frac{h}{3}\Big[f(x_0)+4f(x_1)+f(x_2)\Big] = \frac{(\frac{1}{2})}{3}\Big[3+4\frac{29}{4}+15\Big]=\frac{47}{6} $$ $$
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * style="width:20%; padding:10px; border:2px solid #8888aa" |
 * $$\displaystyle \longrightarrow (2)
 * style= |
 * }.
 * }.

in coclusion, From equation $$\displaystyle (1) ,\ (2)$$,

$$\displaystyle I$$ and $$\displaystyle I_2$$ is same for 3rd polynomial equation


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

$$\displaystyle I=I_2=\frac{47}{6} $$
 * style="width:10%; padding:10px; border:2px solid #8888aa" |
 * style="width:10%; padding:10px; border:2px solid #8888aa" |
 * style= |
 * }.
 * }.

= Problem 16: Runge phenomenon =

Given
Refer Lecture slide 16-1 for problem statement
 * {| style="width:100%" border="0" align="left"

I = \int^{5}_{-5} \frac{1}{1+x^2} dx $$ $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle \longrightarrow (1)


 * }.
 * }.

Find
(1) Find $$\displaystyle I_n$$ using Newton-Cotes for n=1, 2,. . ., 15

(2) Plot $$\displaystyle f, f_n $$ $$\mbox{,  }$$ for n=1, 2, 3, 8, 12

(3) Plot $$\displaystyle I_n$$ vs. $$\displaystyle n$$ and observe that $$\displaystyle I_n$$ does not converge as $$\displaystyle n \rightarrow \infty $$

(4) Observe weight $$\displaystyle W_{i,n}:=\int^{b}_{a} l_{i,n}(x),\ dx$$ are not all positive for n $$ \ge $$ 8. Plot $$\displaystyle l_{i,n}(x)$$ for i=1, 2, ... ,8 & n=8.


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


 * $$\displaystyle
 * $$\displaystyle

$$


 * }.
 * }.

Solution
(1) Find $$\displaystyle I_n$$ using Newton-Cotes for n=1, 2,. . ., 15

Newton-cotes formula |$$\displaystyle P_n$$,


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

P_n = \sum_{i=0}^n l_{i,n}(x)f(x_i)$$
 * $$\displaystyle
 * $$\displaystyle


 * }.
 * }.

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

l_{i,n}(x) = \prod_{j=0,j\ne i}^n \frac{x-x_j}{x_i-x_j}$$
 * $$\displaystyle
 * $$\displaystyle


 * }.
 * }.


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

I_n = \int^{5}_{-5}P_n(x) dx = \sum_{i=0}^n \underbrace{\int^{5}_{-5}l_{i,n}(x) dx}_{W_i} f(x_i) = \sum_{i=0}^n w_i f(x_i) $$ $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle \longrightarrow (2)
 * }.
 * }.

From the lecture slide 13-1,

the Error of integration of Newton-Cotes formula


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

\left | E_n(f(x)) \right | = \left | I-I_n \right | \le \int^{b}_{a} \left | f(x)-P_n(x)\right |  dx \le \frac {\int^{b}_{a}\left | q_{n+1}(x) \right| dx}{(n+1)!}\Big( max f^{(n+1)}(\xi), \xi \in [a,b] \Big)$$ $$
 * $$\displaystyle
 * $$\displaystyle
 * $$\displaystyle \longrightarrow (3)


 * }.
 * }.

From the equation (2),

We know that when the $$\displaystyle f(x) $$ is the polynomial of degree less than $$n$$


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

E_n(f(x)) = 0 $$
 * $$\displaystyle
 * $$\displaystyle


 * }.
 * }.

because, for the polynomial of degree less than $$\displaystyle n$$


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

if\quad f(x) \in P_n, \qquad f^{(n+1)}(x)=0 $$
 * $$\displaystyle
 * $$\displaystyle


 * }.
 * }.

In order to find the $$\displaystyle w_i$$ in equation (2), we can use n equations and solve using matrix method like below,

for the $$\displaystyle f(x)=1,\; x, \; x^2, \; ...... \; x^n$$, we can plug f(x) into equation(2)


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

E_n(f(x))= E_n(1)=E_n(x)=E_n(x^2)=......=E_n(x^n)=0 $$
 * $$\displaystyle
 * $$\displaystyle


 * }.
 * }.


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



\begin{align} &E_n(1)= \int^{5}_{-5}1 dx -[w_1(1)+w_2(1)+w_3(1)+......+w_n(1)] =0\\ &E_n(x)= \int^{5}_{-5}x dx -[w_1(x_0)+w_2(x_2)+w_3(x_3)+......+w_n(x_n)] =0\\ &E_n(x^2)= \int^{5}_{-5}x^2 dx -[w_1(x_0^2)+w_2(x_1^2)+w_3(x_3^2)+......+w_n(x_n)^2] =0\\ &E_n(x^3)= \int^{5}_{-5}x^3 dx -[w_1(x_0^3)+w_2(x_1^3)+w_3(x_3^3)+......+w_n(x_n)^3] =0\\ &.......................\\ &..............................\\ &......................................\\ &E_n(x^n)= \int^{5}_{-5}x^n dx -[w_1(x_0^n)+w_2(x_1^n)+w_3(x_3^n)+......+w_n(x_n)^n] =0 \end{align} $$
 * $$\displaystyle
 * }.
 * }.

Hence, the $$\displaystyle w_i$$ may be determined from


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



\underbrace{ \begin{bmatrix} 1 & 1 & 1 & 1 & ...... & 1\\ x_0 & x_1 & x_2 & x_3 & ...... & x_n\\ x_0^2 & x_1^2 & x_2^2 & x_3^2 & ...... & x_n^2\\ x_0^3 & x_1^3 & x_2^3 & x_3^3 & ...... & x_n^3\\ \vdots & \vdots& \vdots& \vdots& \vdots & \vdots \\ x_0^n & x_1^n & x_2^n & x_3^n & ...... & x_n^n\\ \end{bmatrix} }_{A} \underbrace{ \begin{bmatrix} w_0\\ w_1\\ w_2\\ w_3\\ \vdots\\ w_n\\ \end{bmatrix} }_{W} = \underbrace{ \begin{bmatrix} \int^{b}_{a}1 dx = (b-a)\\ \int^{b}_{a}x dx = (b^2-a^2)/2\\ \int^{b}_{a}x^2dx = (b^3-a^3)/\\ \int^{b}_{a}x^3dx = (b^4-a^4)/4\\ \vdots \\ \int^{b}_{a}x^ndx =(b^n-a^n)/(n+1)\\ \end{bmatrix} = \begin{bmatrix} (b-a)\\ \frac{(b^2-a^2)}{2}\\ \frac{(b^3-a^3)}{3}\\ \frac{(b^4-a^4)}{4}\\ \vdots \\ \frac{(b^{n+1}-a^{n+1})}{n+1} \end{bmatrix} }_{B} $$
 * $$\displaystyle
 * }.
 * }.

$$\displaystyle W_i$$ can be obtained like this,
 * {| style="width:100%" border="0" align="left"

\begin{align} &AW=B\\ &W_i=A^{-1}B \end{align} $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * }.

plug obtained $$\displaystyle W_i$$ and \displaystyle f(x)=\frac{1}{1+x^2} into equation(2), $$\displaystyle I_n$$ can be calculated,
 * {| style="width:100%" border="0" align="left"

I_n = \sum_{i=0}^n w_i f(x_i) $$
 * $$\displaystyle
 * $$\displaystyle
 * }.
 * }.

$$ \begin{array}{|c||c|c|} n & I_n \\ \hline 1&0.384615\\ 2&6.794872\\	3&2.081448\\ 4&2.374005\\ 5&2.307692\\ 6&3.870449\\	7&2.898994\\	8&1.500489\\ 9&2.398618\\ 10&4.673301\\ 11&3.244773\\ 12&-0.31294\\ 13&1.919797\\ 14&7.899545\\ 15&4.155559\\ \hline \end{array} $$

Matlab Code:

(2) Plot $$\displaystyle f, f_n $$ $$\mbox{,  }$$ for n=1, 2, 3, 8, 12

$$\displaystyle f, f_n $$ were plotted by below matlab code.

Matlab Code:

(3) Plot $$\displaystyle I_n$$ vs. $$\displaystyle n$$ and observe that $$\displaystyle I_n$$ does not converge as $$\displaystyle n \rightarrow \infty $$



(4) Observe weight $$\displaystyle W_{i,n}:=\int^{b}_{a} l_{i,n}(x),\ dx$$ are not all positive for n $$ \ge $$ 8. Plot $$\displaystyle l_{i,n}(x)$$ for i=1, 2, ... ,8 & n=8.

This is the table that show the $$\displaystyle W_{i,n}:=\int^{b}_{a} l_{i,n}(x),\ dx$$ for n=8


 * Plot $$\displaystyle l_{i,n}(x)$$ for i=1, 2, ... ,8 & n=8.
 * Plot $$\displaystyle l_{i,n}(x)$$ for i=1, 2, ... ,8 & n=8.


 * Matlab Code: