User:Egm6341.s10.team3.ynahn81/HW4-8(c)

(c) Compare Comp Trap, Comp Simp, Romberg for $$\displaystyle f(\theta) = \frac{1-sin^{2}(\frac{\pi}{12})}{1-sin(\frac{\pi}{12})cos\theta}$$
Ref: Lecture notes [[media:Egm6341.s10.mtg25.djvu|p.25-2]]

Problem Statement
For $$\displaystyle I = \int_{0}^{2\pi}\frac{1-sin^{2}(\frac{\pi}{12})}{1-sin(\frac{\pi}{12})cos\theta}d\theta$$, compare $$\displaystyle I_{n}$$ to $$\displaystyle O(10^{-10})$$ and computational time for - Composite Trap rule, - Composite Simpson's rule and - Romberg table.

Solution
In what follows, we are using the follow matlab function as a given function $$\displaystyle f(\theta) = \frac{1-sin^{2}(\frac{\pi}{12})}{1-sin(\frac{\pi}{12})cos\theta}$$.

i) Composit Trap. rule

Recall the formula of comp. Trap. rule.


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

T_{0}(n) = h\left[\frac{1}{2}f(x_{0}) + f(x_{1}) + f(x_{2}) + \cdot\cdot\cdot + f(x_{n-1}) + \frac{1}{2}f(x_{n}) \right] $$ where $$\displaystyle h = \frac{b-a}{n} $$ and $$\displaystyle x_{0} = a,\quad x_{n} = b $$.
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Also, realize the fact that if we want to double the number of sub-intervals, we don't have to calculate the functions at even point, i.e. $$\displaystyle f(x_{0}), f(x_{2}), \cdot\cdot\cdot, f(x_{n-2}), f(x_{n}) $$ because we already calculated the functions at these even point at the previous step. Therefore, all we have to do is evaluating the functions only at the odd points and combine them with the calculated values at the previous step. As a result, we can use the following inductive formula.
 * {| style="width:100%" border="0" align="left"

T_{0}(2n) = T_{0}(n) + h\sum_{i=1}^{n}f(x_{2i-1}). $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Matlab code

Result

$$\displaystyle I_{n,trap} = T_{0}(n) = 6.069090959564874, \quad n = 2^{4} = 16 $$ $$\displaystyle E_{n,trap} = I - I_{n,trap} = 1.039168751049147e-13 $$

ii) Composite Simpson's rule

Recall the formula of comp. Simpson's rule.


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

T_{1}(n) = \frac{h}{3}\left[f(x_{0}) + 4f(x_{1}) + 2f(x_{2}) + 4f(x_{3}) + \cdot\cdot\cdot + 4f(x_{n-3}) + 2f(x_{n-2}) + 4f(x_{n-1}) + f(x_{n}) \right] $$ where $$\displaystyle h = \frac{b-a}{n} $$ and $$\displaystyle x_{0} = a,\quad x_{n} = b $$.
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Also, realize the fact that if we want to double the number of sub-intervals, we don't have to calculate the functions at even point, i.e. $$\displaystyle f(x_{0}), f(x_{2}), \cdot\cdot\cdot, f(x_{n-2}), f(x_{n}) $$ because we already calculated the functions at these even point at the previous step. Therefore, all we have to do is evaluating the functions only at the odd points and combine them with the calculated values at the previous step. As a result, we can use the following inductive formula.
 * {| style="width:100%" border="0" align="left"

T_{1}(2n) = T_{1}(n) + h\sum_{i=1}^{n}4f(x_{2i-1}) - \underbrace{h\sum_{i=1}^{\frac{n}{4}}2f(x_{4i-2})}_{Already\,\, obtained\,\, when\,\, T_{1}(n) \,\, was\,\, calculated.}. $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Matlab code

Result

$$\displaystyle I_{n,Simpson} = T_{1}(n) = 6.069090959564741, \quad n = 2^{5} = 32 $$ $$\displaystyle E_{n,Simpson} = I - I_{n,Simpson} = 2.930988785010413e-14 $$

iii) Romberg Table

Since the first column of Romberg table corresponds to Trap. rule, we can use the below formulas in order to obtain the first column of Romberg talbe.


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

T_{0}(n) = h\left[\frac{1}{2}f(x_{0}) + f(x_{1}) + f(x_{2}) + \cdot\cdot\cdot + f(x_{n-1}) + \frac{1}{2}f(x_{n}) \right] $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }


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

T_{0}(2n) = T_{0}(n) + h\sum_{i=1}^{n}f(x_{2i-1}). $$ where $$\displaystyle h = \frac{b-a}{n} $$ and $$\displaystyle x_{0} = a,\quad x_{n} = b $$.
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Also, we gradually set up the Romberg table in parallel, as soon as the result for $$\displaystyle T_0(2^j) $$ becomes available using the below formula.


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

T_{k}(2n) = \frac{2^{2k}T_{k-1}(2n) - T_{k-1}(n)}{2^{2k}-1}. $$
 * $$\displaystyle
 * $$\displaystyle
 * }
 * }

Matlab code

Result

$$\displaystyle I_{n,Romberg} = T_{0}(n) = 6.069090959564874, \quad n = 2^{4} = 16 $$ $$\displaystyle E_{n,Romberg} = I - I_{n,Romberg} = 1.039168751049147e-13 $$

iv) Compare computational time

We can compare computational times of above three methods using tic/toc commands in Matlab.

Matlab code

As above table shown, Simpson's rule has the shortest computational time. However, we should note that the computational time does not really different from each other for three method. Especially, trap rule has a shorter computational time than Romberg Table (in general, computational time of Romberg Table is much shorter than that of Trap. rule.). Since the given function $$\displaystyle f(\theta) = \frac{1-sin^{2}(\frac{\pi}{12})}{1-sin(\frac{\pi}{12})cos\theta}$$ is a periodic function, the numerical integration of this function is converged very fast by Trap. rule. Therefore, '''Trap. rule is an efficient way to numerically integrate the given function $$\displaystyle f(\theta) = \frac{1-sin^{2}(\frac{\pi}{12})}{1-sin(\frac{\pi}{12})cos\theta}$$ because this function is a periodic function.'''