University of Florida/Egm6341/s10.Team2/HW3

Problem Statement
-Pg. 17-1


 * Derive the tighter error bound of the simple Simpson's rule on pg.14-2 for the following two cases by changing the function:

$$ G\left ( t \right ):= e\left ( t \right )-t^{4} \; e\left ( 1 \right )$$ $$ G\left ( t \right ):= e\left ( t \right )-t^{6} \; e\left ( 1 \right )$$

Point out where the proof breaks down.


 * $$ G\left ( t \right ):= e\left ( t \right )-t^{5} \; e\left ( 1 \right )$$ Find if $$ G^{\left (3 \right )}\left ( 0 \right )= 0$$ and follow the same steps as in the proof to see what happens.

Derivation of Tighter Error Bound of Simple Simpson's Rule
$$\begin{align} E_{2} & = \frac{(b-a)^5}{2880}f^{(4)}(\xi )\\ & = \frac{h^5}{2880}f^{(4)}(\xi )\\ \end{align}$$

$$ \mbox{where, }\quad h=\frac{(b-a)}{2}$$

Shift origin of x axis to point at $$x_1\!$$

$$x(t)=x_1+ht\! \qquad \mbox{where, } t\epsilon [-1,1]\!$$

$$\begin{cases} t=0 &, x= x(0) = x_1\\ t=+1 &, x= x(1) = x_1+h = x_2\\ t=-1 &, x= x(-1)= x_0\\ \end{cases}$$

$$E_2= I - I_2 = \frac{(b-a)}{2} e(1) = h\;e(1)$$

Therefore, $$\;\;-\frac{h^5}{90}f^{(4)}(\xi)= h\;e(1)$$

$$e(1)=\frac{h^4}{90}f^{(4)}(\xi)$$

$$e(t)=\int_{-t}^{+t}F(t)dt-\frac{t}{3}[F(-t)+4F(0)+F(t)] \qquad \mbox{where, } f(x(t))=F(t)$$

By Definition: $$\quad G^{(3)}(\xi_3)=0$$


 * $$ G\left ( t \right ):= e\left ( t \right )-t^{4} \; e\left ( 1 \right )$$

$$G^{(3)}(t)=e^{(3)}(t)-\frac{t^7}{210}\;e(1)$$

$$G^{(3)}(\xi_3)= 0 = e^{(3)}(\xi _3)-\frac{{\xi_{3}}^{7}}{210}\;e(1)$$

$$\begin{align} 0      & = e^{(3)}(\xi _3)-\frac{{\xi_{3}}^{7}}{210}\;e(1)\\ e(1)   & =e^{(3)}(\xi_3)\;\frac{210}{{\xi_{3}}^{7}}\\ \end{align}$$

Proof breaks down here. Neglects to properly integrate $$t^5\;e(1)$$


 * $$ G\left ( t \right ):= e\left ( t \right )-t^{6} \; e\left ( 1 \right )$$

$$G^{(3)}(t)=e^{(3)}(t)-\frac{t^9}{504}\;e(1)$$

$$G^{(3)}(\xi_3)= 0 = e^{(3)}(\xi _3)-\frac{{\xi_{3}}^{9}}{504}\;e(1)$$

$$\begin{align} 0      & = e^{(3)}(\xi _3)-\frac{{\xi_{3}}^{9}}{504}\;e(1)\\ e(1)   & =e^{(3)}(\xi_3)\;\frac{504}{{\xi_{3}}^{9}}\\ \end{align}$$

Proof breaks down here. Neglects to properly integrate $$t^5\;e(1)$$

Solution of Tighter Error Bound Function at t=0
$$ G\left ( t \right ):= e\left ( t \right )-t^{5} \; e\left ( 1 \right )$$

$$G^{(3)}(t)=e^{(3)}(t)-\frac{t^8}{336}\;e(1)$$

$$G^{(3)}(0)=e^{(3)}(0)-\frac{0^8}{336}\;e(1)$$

$$ G^{(3)}(0)= -\frac{2(0)^2}{3}\;[F^{(4)}(0)+90\;e(1)]$$

$$G^{(3)}(0) =0\!$$

Solution for problem 1:Egm6341.s10.team2.patodon 10:05, 17 February 2010 (UTC)

Proofread problem 1:

Problem Statement
Pg. 17-2 Show that the error of the Composite Simpson's Rule is:

$$ \left| E_{n}^{2} \right| \leq \frac{(b-a)^{5}}{2880 \; n^{4}}M_{4} = \frac{(b-a)h^{4}}{2880}M_{4} $$

where $$ \; M_{4} $$ is the maximum of  $$ \quad \left| F^{4}(\zeta) \right| \quad and \quad \zeta \epsilon \left[a,b \right] $$

Problem Solution
The error for the Simpson's Composite rule, can be expressed as follows:

$$ \; E_{n}^{2} = I - I_{n} $$

This error can be rewritten as follows:

$$ E_{n}^{2} =  \int_{a}^{b} f(x)dx - \frac{h}{3}\left[f_{0} + 4f_{1} + 2f_{2} +4f_{3} + 2 f_{4} ... + 2f_{n-2} + 4f_{n-1} + f_{n} \right] $$

This can be expressed by using the Composite Simpson's rule and dividing the interval 'a' to 'b' into smaller intervals, the error can be expressed as follows:

$$ E_{n}^{2} = \int_{x_{i}-1}^{x_{i}} f(x)dx - \frac{h}{3} \left[f(x_{i-1}) + f(x_{i}) \right] $$

In this manner the error for simple Simpson's rule can be used to describe the error as follows:

$$ \left|E_{n}^{2} \right| \leq \left| \frac{h^{5}}{90} \underbrace{\sum_{i=1}^{n} max F^{(4)}(\zeta)}_{\bar{M_{4}}} \right| $$

$$ \bar{M_{4}} \leq n \cdot M_{4} $$

Using this expression the error can then be expressed as follows:

$$ \left|E_{n}^{2} \right| \leq \left| \frac{(b-a)^{5}}{2880 \; n^{5}}n \; M_{4} \right| = \frac{(b-a)^{5}}{2880\; n^{4}}M_{4} $$

$$ \left|E_{n}^{2} \right| \leq \left| \frac{(b-a)(b-a)^{4}}{2880 \; n^{4}}M_{4} \right| = \frac{(b-a)h^{4}}{2880}M_{4} $$

Solution for problem 2: Guillermo Varela

Proofread problem 2: Srikanth Madala 05:22, 18 February 2010 (UTC)

Problem Statement
Pg. 17-2

See also HW1-problem 8HW1-problem8
 * 1) Using the error estimates of Taylor Series ,Composite Trapezoidal rule and Composite Simpson's rule estimate $$\boldsymbol{n}$$ such that $$\boldsymbol{E_n=I-I_n=O(10^{-6})}$$ and compare against numerical results of HW1.
 * 2) Numerically find the power of the step size h in the error i.e Plot the error vs h on semilog and fit a straight line and measure the slope.

Note: This problem is a logical continuation of problem 8 in HW1. In HW1 it was attempted by iteration to find the value of n. In this problem it is done analytically.All numerical data can be found on the HW1 page linked above and also the corresponding MATLAB codes.

Solution: Taylor series
Error is defined as
 * {| style="width:100%" border="0"

$$  \displaystyle E_n=I-I_n =\int_{a}^{b}[f(x)-f_n(x)]dx $$
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

For the Taylor series, from the discussion on p6-4, we know that the error is nothing but the remainder of the Taylor series integrated over the given interval i.e
 * {| style="width:100%" border="0"

$$  \displaystyle \int_{a}^{b}R_{n+1}(x)= \int_{a}^{b}\left [\frac{(x-x_0)^{n+1}}{(n+1)!}f^{(n+1)}(\xi ) \right ] $$ where $$\xi \varepsilon [a,b] $$    (1 p2-3)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

For the given function the error is
 * {| style="width:100%" border="0"

$$  \displaystyle E_n=\int_{0}^{1}\underbrace{\frac{x^n}{(n+1)!}}_{w(x)}\underbrace{e^{\xi (x)}}_{g(x)}dx $$    (2)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

Using the Integral Mean Value Theorem, we have
 * {| style="width:100%" border="0"

$$  \displaystyle E_n=e^{\xi (x=\alpha)} \int_{0}^{1}\frac{x^n}{(n+1)!}dx $$    (3)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

integrating we get,
 * {| style="width:100%" border="0"

$$  \displaystyle E_n=e^{\xi (x=\alpha)} \frac{1}{(n+1)!(n+1)} $$    (4)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

This function has a minimum when $$\; \alpha=0$$ and maximum when $$\; \alpha=1$$ thus we have
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{1}{(n+1)!(n+1)}\leq E_n\leq \frac{e}{(n+1)!(n+1)} $$    (5) Setting the upper bound of the error to $$\; 10^{-6}$$ we have
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }
 * {| style="width:100%" border="0"

$$  \displaystyle \ E_n\leq \frac{e}{(n+1)!(n+1)}=10^{-6} $$    (6) Solving this equation for $$\; n$$ we get
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }


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

$$  \displaystyle (n+1)!(n+1)=e*10^6 \Rightarrow n\approx 7.92 = 8 $$     (A) 
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * 
 * }

In the numerical analysis in problem 8 of HW1 it was seen that the Error for n = 8 was -3.848E-07

Solution:Composite Trapezoidal rule
From the discussion on page 16-3 we have


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

$$  \displaystyle \left | E_n^1 \right |\leq \frac{(b-a)}{12n^2}M_2 $$    (1)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

where $$\; M_2 = max \left | f^{(2)} (\zeta )\right |$$ for $$\; \zeta \varepsilon [a,b]$$

For the given function $$\; f(x)= \frac{e^x - 1}{x}$$, we have


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

$$  \displaystyle f^{(2)}(x)= \frac{e^x [x^2-2x+2]-2}{x^3} $$     (2)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

Evaluating over the given interval it is seen that the function $$\; f^{(2)}(x)$$ has maximum value at $$\; x=1$$


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

$$  \displaystyle M_2= f^{(2)}(x=1)=\frac{e[1+2-2]-2}{1}=0.718282 $$    (3)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

Setting the error to the $$\; 10^{-6}$$ and solving for $$\; n$$ we get


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

$$  \displaystyle \frac{(1-0)^3}{12n^2}(0.718282)=10^{-6}\Rightarrow n\approx 244.657= 245 $$     (B) 
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * 
 * }

In the numerical analysis the value of n =245 was found to be 5.872E-07

Solution :Composite Simpson's Rule
We have from p 17-2, the error estimate of the Composite Simpson's rule as
 * {| style="width:100%" border="0"

$$  \displaystyle \left | E_n^2 \right |\leq \frac{(b-a)^5}{2880n^4}M_4 $$    (1)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

where $$M_4 = max \left | f^{(4)} (\zeta )\right |$$ for $$\zeta \varepsilon [a,b]$$

For the given function $$f(x)= \frac{e^x - 1}{x}$$, we have


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

$$  \displaystyle f^{(4)}(x)= \frac{e^x [x^4-4x^3+12x^2-24x]-24}{x^5} $$     (2)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

Evaluating over the given interval it is seen that the function $$\; f^{(2)}(x)$$ has maximum value at $$\; x=1$$
 * {| style="width:100%" border="0"

$$  \displaystyle M_4= f^{(4)}(x=1)=\frac{e[1-4+12-24+24]-24}{1}=9(e)-24=0.464536
 * style="width:95%" |
 * style="width:95%" |

$$    (3)
 * <p style="text-align:right">
 * }

Setting the error to the $$\; 10^{-6}$$ and solving for $$\; n$$ we get


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

$$  \displaystyle \frac{(1-0)^5}{2880n^4}(0.464536)=10^{-6}\Rightarrow n\approx 3.35735= 4 $$     (c) <div style="width: 80%; margin-left: auto; margin-right: auto; padding: 4px; border: 2px solid #FF0000; background-color: #FFDDDD; text-align: left;">
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

In the numerical analysis the error for n= 4 was 6.717E-06

Part 2:Numerical determination of power of h
In order to verify the power of $$\boldsymbol {h}$$ in the error, data from Problem 8 of HW1 is used.

In the case of a Semilog plot [log(y) vs x], a straight line has an equation of the form

$$\boldsymbol{log(Y)=a X+b}$$

such that $$\boldsymbol{Y= e^{(aX+b)}=e^b (e^{aX})=ke^{aX}}$$.From the above equation it is seen that if the plot on a semilog graph is a straight line then the relationship between the two variables is exponential.

A log-log plot(log(y) vs log(x)) is a similar plot, for which the equation is of the form

$$\boldsymbol{log(Y)=a [log(X)]+b}$$

such that, $$\boldsymbol{Y= kx^a;k=e^b}$$.It is readily seen that the slope of the line in the log-log graph is the power of the x variable. Further reference on the theory of semilog and log plots and methods to fit curves is given below:
 * 1) Semilog plot [w:Semi-log_graph]
 * 2) Log-Log plot [w:Log-log_graph]
 * 3) Methods of fitting curves in MATLAB for power and exponential series []

This discussion is used in the interpretation of the graphs given below.

Composite trapezoidal Rule
Given below is the data from the numerical evaluation of the given function using Composite Trapezoidal Rule.[Egm6341.s10.Team2/HW1] First we plot a semilog graph for the data. The graph is shown below:



The graph is not a straight line which implies that the relationship between $$\boldsymbol{h}$$ and $$\boldsymbol{log(error)}$$ is not exponential. Hence we plot a log-log graph as below:



This is seen to be linear. A straight line if fitted to the data the equation of which is given above. From the discussion above, we see that the slope of line is 2.1447 which is very close to the analytical value of 2.

Composite Simpson's Rule
Given below is the data obtained from HW1 for the Composite Simpsons Rule.[Egm6341.s10.Team2/HW1]

Plotting a semilog graph of $$\boldsymbol{log(error)}$$ against $$\boldsymbol{h}$$ we see that it is non-linear as in the case of the Composite Trapezoidal Rule.



Thus, plotting the log-log graph as below,



we see that the slope of the line is 4.0881 which is very close to the analytically determined value of 4.

Taylor Series
Given below is the data for the Taylor series method.$$\boldsymbol{h}$$ can be defined for the taylor series as $$(b-a)/n$$ but since $$\boldsymbol{h}$$ is not used explicitly in the Taylor series method,it might be of interest to observe how error depends on $$\boldsymbol{n}$$,the number of terms of the series. Given below are the plots.Semilog plots and log-log plots can be made of error vs h much like the ones given above, but no linearity is shown in both.



From this plot it is clearly seen that minimizing the error does not require an infinite series.The error curve is asymptotic and it can be seen that 8 terms gives a very low error.



Solution for problem 3: Niki Nachapa

Proofread problem 3: Guillermo Varela

Problem Statement
Pg. 17-3 and 18-1

a] Replicate the table 5.1, page 255, Atkinson text book using the composite Trapezoidal rule.

$$I=\int_{0}^{\pi} e^{x} \cos \left ( x \right )\cdot dx$$

b] Replicate the table 5.3, page 258, Atkinson text book using the composite Simpson's rule.

$$I=\int_{0}^{\pi} e^{x} \cos \left ( x \right )\cdot dx$$

c] Replicate the table 5.4, page 261, Atkinson text book

$$I=\int_{0}^{1} x^{3} \sqrt{x} \cdot dx$$

d] Replicate the table 5.5, page 262, Atkinson text book

$$I=\int_{0}^{5} \frac{1}{1+\left ( x-\pi \right )^{2}} \cdot dx$$

e] Replicate the table 5.6, page 262, Atkinson text book

$$I=\int_{0}^{1} \sqrt{x} \cdot dx$$

f] Replicate the table 5.7, page 263, Atkinson text book

$$I=\int_{0}^{2\pi} e^{\cos \left ( x \right )} \cdot dx$$

Problem Solution
 Solution for Part (a):  Here our given Integral is: $$I=\int_{0}^{\pi} e^{x} \cos \left ( x \right )\cdot dx$$ where $$f\left ( x \right )=e^{x} \cos \left ( x \right )$$ is the integrand and $$\;x\in \left [ 0,\pi \right ]$$. We use composite Trapezoidal rule to obtain the following tabulated results:

 MATLAB Code: 

 Solution for Part (b):  Here our given Integral is: $$I=\int_{0}^{\pi} e^{x} \cos \left ( x \right )\cdot dx$$ where $$f\left ( x \right )=e^{x} \cos \left ( x \right )$$ is the integrand and $$\;x\in \left [ 0,\pi \right ]$$. We use composite Simpson's rule to obtain the following tabulated results:

 MATLAB Code: 

 Solution for Part (c):  Here our given Integral is: $$I=\int_{0}^{1} x^{3} \sqrt{x} \cdot dx$$ where $$f\left ( x \right )=x^{3} \sqrt{x}$$ is the integrand and $$\;x\in \left [ 0,1 \right ]$$. We use both composite Trapezoidal and Simpson's rules to obtain the following tabulated results:

 MATLAB Code: 

 Solution for Part (d):  Here our given Integral is: $$I=\int_{0}^{5} \frac{1}{1+\left ( x-\pi \right )^{2}} \cdot dx$$ where $$f\left ( x \right )=\frac{1}{1+\left ( x-\pi \right )^{2}}$$ is the integrand and $$\;x\in \left [ 0,5 \right ]$$. We use both composite Trapezoidal and Simpson's rules to obtain the following tabulated results:

 MATLAB Code: 

 Solution for Part (e):  Here our given Integral is: $$I=\int_{0}^{1} \sqrt{x} \cdot dx$$ where $$f\left ( x \right )=\sqrt{x}$$ is the integrand and $$\;x\in \left [ 0,1 \right ]$$. We use both composite Trapezoidal and Simpson's rules to obtain the following tabulated results:

 MATLAB Code: 

 Solution for Part (f):  Here our given Integral is: $$I=\int_{0}^{2\pi} e^{\cos \left ( x \right )} \cdot dx$$ where $$f\left ( x \right )=e^{\cos \left ( x \right )}$$ is the integrand and $$\;x\in \left [ 0,2\pi \right ]$$. We use both composite Trapezoidal and Simpson's rules to obtain the following tabulated results:

 MATLAB Code: 

Solution for problem 4: Srikanth Madala 05:22, 18 February 2010 (UTC)

Proofread problem 4:

Problem Statement
Pg. 19-2 (cont'd on Pg. 6-5)

a) Modify the code to make the computation of T(2^j) efficient;

b) Add Romber Table and compare with previous results.

Problem Solution
a) Modified code

function [I1,I2,T1] = trapzoid(f,a,b,n) h=(b-a)/n; x=linspace(a,b,n+1); fx=feval(f,x); I1=h*(fx(1)/2+sum(fx(2:1:n))+fx(n+1)/2); n=2*n; h=(b-a)/n; x=linspace(a,b,n+1); fx=feval(f,x); I2=I1/2+h*sum(fx(2:2:n)); T1=(2^2*I2-I1)/(2^2-1);

%call function %I1=To(4),I2=To(8), and T1=T1(4) format long [I1,I2,T1]=trapzoid(@(x) (exp(x)-1)./x,eps,1,4)

b) Romberg Table



Convergence is achieved much faster than the composite trapezoidal rule.

Solution for problem 5:

Proofread problem 5:

Contributing Authors
Problem 1 Solution --Patrick O'Donoughue 10:07, 17 February 2010 (UTC)

Problem 2 solution --Guillermo Varela 19:16, 17 February 2010 (UTC)

Problem 3 solution --Niki Nachappa 15:31, 17 February 2010 (UTC)

Problem 4 solution and proof-read problem 2 -- Srikanth Madala 05:22, 18 February 2010 (UTC)

Problem 5 solution -- PENGXIANG JIANG 19:00, 17 February 2010 (UTC)